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ABSTRACT 


In  this  thesis  an  investigation  is  made  of  the 
properties  of  long  aeromagnetic  profiles  in  order  to  deter¬ 
mine  some  aspects  of  the  magnetic  character  of  the  crust  of 
the  Earth.  In  particular,  an  interpretation  is  made  of  those 
broad  anomalies  whose  widths  lie  between  100  and  150  miles. 
These  broad  anomalies  are  extracted  from  the  profile  data 
using  numerical  filtering.  The  frequency  characteristics  of 
the  filters  to  be  applied  are  determined  by  performing  a 
power  spectral  analysis  on  the  data.  The  spectral  analyses 
also  provide  some  insight  as  to  the  magnetic  character  of 
the  crust. 

The  many  laborious  calculations  involved  in  such 
an  investigation  have  been  evaluated  with  the  aid  of  the 
IBM  1620  digital  computer  in  the  Computing  Centre,  University 
of  Alberta. 

The  broad  anomalies  revealed  by  numerical  filtering 
are  interpreted  as  being  due  to  anomalous  Curie  point  geotherm 
structure  and  estimates  of  the  depth  to  the  geotherm  and  of 
susceptibility  contrast  at  these  depths  are  made  for  both  the 
continental  and  oceanic  crusts. 
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The  results  suggest  that: 

(l)  the  Curie  point  depth  for  the  continental  crust  is 
slightly  greater  than  that  for  the  oceanic  crust, 

(ii)  the  geotherm  depth  is  very  close  to  that  of  the 
Mohorovicic  discontinuity  for  the  continental  crust 
and  much  greater  than  that  of  the  discontinuity  for 
the  oceanic  crust, 

(iii)  the  susceptibilities  of  rocks  at  the  geotherm  depth 
are  greater  for  the  continental  crust. 
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I.  INTRODUCTION 


1 . 1  The  Problem 

In  many  recent  geophysical  publications  (Dobrin 
and  Ward,  1962;  Goldstein,  1962;  Dean,  1958)  it  has  been 
suggested  that  statistical  analysis  and  the  application  of 
filtering  to  geophysical  data  might  reveal  a  new  technique 
in  the  interpretation  of  many  problems,  especially  in  problems 
of  potential  field  signals.  Although  the  mathematics  for 
such  operations  have  been  known  for  two  or  three  decades, 
the  laborious  calculations  required  for  such  interpretation 
were  never  thought  to  be  practical  and  too,  there  were  those, 
as  there  still  are  today,  who  felt  that  potential  field 
problems  do  not  lend  themselves  to  such  an  analysis.  With 
the  advent  in  recent  years  of  high-speed  analog  and  digital 
computers  much  interest  and  research  has  been  initiated 
along  these  lines  in  hope  that  some  light  might  be  shed  upon 
problems  whose  solutions  have  escaped  other  interpretational 
procedure  or  have  lacked  corroboration  by  a  more  detailed 
type  of  analysis. 

One  such  problem  has  been  that  of  the  Curie  point 
geotherm,  the  depth  to  which  is  very  important  in  the  inter¬ 
pretation  of  deep  magnetic  sources.  Such  information  could 
also  give  the  interpreter  an  insight  as  to  the  rock  type 
which  occurs  at  the  geotherm  depth. 
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In  the  past,  the  Curie  point  geotherm  has  always 
been  interpreted  as  the  lower  boundary  of  magnetic  sources 
which  give  rise  to  the  broader  features  of  magnetic  maps  or 
profiles.  In  this  thesis  it  is  shown  that  certain  broad 
anomalies  of  an  aeromagnetic  profile  (with  dipole  field 
removed)  as  revealed  by  spectral  analysis  and  numerical 
filtering  could  be  due  to  upwarping  or  downwarping  of  the 
geotherm.  The  analysis  is  performed  on  data  from  both  contin¬ 
ental  and  oceanic  aeromagnetic  profiles  and  it  is  hoped  that 
some  insight  will  be  gained  into  the  difference  between 
continental  and  oceanic  magnetic  crustal  character  as  well 
as  the  nature  of  the  Curie  point  geotherm. 

1 . 2  Terrestrial  Magnetism 

The  existence  of  a  magnetic  field  associated  with 
the  earth  was  discovered  as  early  as  the  eleventh  century 
when  mariners  began  using  the  magnetic  compass  in  navigation. 
However,  it  was  not  until  l600,  that  the  true  nature  of  the 
magnetic  field,  as  observed  at  the  earth's  surface,  was 
revealed  by  the  experimental  work  of  William  Gilbert.  His 
publication  "De  Magnete"  pointed  out  that  the  earth  itself 
acts  to  a  certain  extent  as  though  it  were  a  great  spherical 
magnet.  Measurements  of  the  magnetic  field  and  its  varia¬ 
tions  with  time  have  been  made  since  1600  and  "Terrestrial 
Magnetism"  .may  therefore  be  described  as  one  of  the  oldest 
fields  in  geophysics. 
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Investigation  of  the  earth’s  magnetic  field 
utilizing  potential  theory  has  led  certain  observers  to 
believe  that  it  originates  from  sources  inside  the  earthy 
outside  the  earth.,  and  on  the  surface  of  the  earth  itself. 

By  applying  spherical  harmonics  to  observations 
of  the  field  made  at  the  earth’s  surface  it  has  been  found 
that  9^  Pei5  cent  of  the  field  originates  from  sources  inside 
the  earth.  Furthermore  the  observed  internal  field  can  be 
described  as  the  effect  of  a  number  of  fictitious  magnetic 
dipoles,  each  with  a  different  orientation,  located  at  the 
centre  of  the  earth.  The  first  of  these  imaginary  dipoles 
is  located  with  its  axis  along  the  magnetic  axis  of  the 
earth  whose  poles  are  displaced  by  approximately  18°  of  lati¬ 
tude.  Because  they  are  not  diametrically  opposite  one 
another  the  poles  are  joined  by  a  line  which  passes  some  750 
miles  from  the  earth's  centre.  Because  of  the  Irregularities 
in  the  earth’s  magnetic  field,  as  observed  at  the  surface, 
an  indefinite  number  of  weaker  dipoles  with  different  orienta¬ 
tions  are  added;  their  various  strengths  and  orientations 
are  determined  from  global  measurements  of  the  field  such 
that  the  best  possible  fit  to  all  the  data  is  obtained. 
However,  if  one  assumes  the  earth  to  be  uniformly  magnetized 
and  that  only  the  axially  oriented  dipole  is  present  then 
the  value  of  polarization  required  to  account  for  the  earth’s 
field  is  many  orders  of  magnitude  higher  than  that  of 
igneous  rocks  associated  with  the  earth’s  crust.  This, 
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coupled  with  the  fact  that  as  the  temperature  of  a  permanently 
magnetized  material  Is  increased  it  loses  its  ability  to 
become  permanentaly  magnetized,  make  this  theory  of  imaginary 
dipoles  largely  untenable  since  the  Curie  points1  of  all 
magnetic  materials  are  reached  at  very  shallow  depths  below 
the  earth’s  surface. 

Experimentation  on  the  effect  of  increased  pressure 
on  the  Curie  point  of  some  materials  shows  that  the  Curie 
temperature  can  be  increased  but  in  general  it  is  believed 
that  the  change  would  not  be  sufficient  to  make  permanent 
magnetization  possible  within  the  depths  of  the  earth. 

The  most  widely  held  theory  as  to  the  origin  of 
the  internal  field  is  that  it  is  due  to  motions  within  the 
liquid  core.  Two  possible  causes  of  motion  within  the  core 
have  been  proposed.  The  first  ascribes  it  to  the  variation 
in  the  angular  velocity  of  the  earth  due  to  tidal  friction 
and  the  second  is  derived  from  a  temperature  gradient  within 
the  core  sufficiently  large  enough  to  produce  thermal  con¬ 
vection.  Because  the  first  effect  is  extremely  small  it 
has  been  ruled  out  as  a  possible  cause  of  the  required 
motion.  Upon  acceptance  of  the  second  cause  then  speculation 
was  made  as  to  the  mechanism  under  which  motions  in  the  core 
might  set  up  electromotive  forces. 


Curie  point  -  the  temperature  at  which  a  substance  upon 

being  heated  loses  its  ferromagnetic  proper¬ 
ties  and  begins  to  behave  as  a  paramagnetic 
substance  if  placed  in  a  magnetic  field. 
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Again  two  ways  have  been  put  forward.  The  first 
is  that  the  motions  set  up  temperature  distributions  within 
the  core  giving  rise  to  thermoelectric  currents.  According 
to  Runcorn  (1956)  potential  differences  across  a  contact 
between  two  spherical  surfaces  of  dissimilar  materials  or 
similar  materials  under  different  physical  conditions  will 
excite  the  electric  mode.  Runcorn  further  states  that  a 
potential  difference  of  1  volt  is  sufficient  to  produce 
sufficient  electric  mode  magnetic  field.  Quantitative 
results  are  not  available,  however,  for  this  theory  at 
present . 

The  second  means  by  which  motion  might  set  up 
electromotive  forces  in  the  core  is  known  as  the  dynamo 
theory.  This  theory  suggests  that  the  field  is  produced 
and  maintained  by  an  induction  mechanism  with  the  kinetic 
energy  of  the  fluid  motions  in  the  core  being  converted  to 
magnetic  energy.  The  crux  of  the  objection  to  this  hypothesis 
is  whether  a  sphere,  possessing  asymmetric  fluid  motion,  but 
unlike  a  real  dynamo  not  possessing  asymmetry  of  structure 
can  in  fact  be  a  self-exciting  dynamo  independent  of  a 
sustaining  field  from  an  external  source.  Bullard  and 
Gellman  (195^)^  using  one  particular  set  of  motions,  have 
devised  a  solution;  the  convergence  of  this  solution  how¬ 
ever  is  debatable.  Nevertheless,  there  seems  to  be  a  general 
consensus  of  opinion  that  motions  in  the  fluid  core  are 
primarily  responsible  for  the  internal  magnetic  field  of  the 


earth . 
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Only  a  small  percentage  of  the  earth’s  field  need 
be  accounted  for  by  an  external  field.  Most  of  the  theories 
proposed  to  explain  the  source  of  the  external  field  Involve 
the  inductive  effect  of  circulating  electric  currents  in  the 
ionosphere . 

The  component  of  the  earth's  field  purported  to  lie 
on  the  surface  of  the  earth  and  sometimes  called  the  non¬ 
potential  field  is  supposedly  due  to  electric  currents  flowing 
from  the  earth  to  the  atmosphere  and  back  again  across  the 
earth’s  surface.  However ,  many  feel  that  the  residue  of 
field  for  which  this  theory  attempts  to  account  is  actually 
due  to  errors  in  observation  and  not  to  any  source  at  all. 

The  earth’s  magnetic  field  has  associated  with  it 
variations  with  time  in  both  strength  and  direction  .  These 
variations  can  be  categorized  under  two  main  types:  secular 
variations  or  changes  in  the  earth's  field  progressing  over 
decades  or  centuries  and  diurnal  variations  or  daily  changes 
which  are  smaller  but  more  rapid  in  oscillation.  A  very 
comprehensive  treatment  of  these  variations  can  be  found  in 
Chapman  and  Bartels  treatise  on  geomagnetism. 

1 . 3  Magnetization  in  the  Crust 

Having  described  the  terrestrial  magnetic  field  and 
its  possible  origin  an  examination  will  now  be  made  of  the 
effect  of  this  field  on  the  magnetically  susceptible  materials 
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In  the  earth's  crust.  With  respect  to  local  areas,,  the 
earth's  magnetic  field  can  be  said  to  be  uniform  over  areas 
of  magnetically  homogeneous  composition  in  the  crust,  but 
is  measurably  distorted  in  regions  of  inhomogeneous  composi¬ 
tion.  Most  of  these  magnetic  anomalies  are  due  to  igneous 
rocks,  sedimentary  rocks  which  contain  iron  oxides  and  to 
some  extent  concentrations  of  magnetic  iron  ores.  The  strength 
of  an  anomaly  generally  depends  on  the  susceptibilities'*'  of 
the  materials  in  the  subsurface  structure  and  also  on  the 
extent  of  the  structure . 

Generally,  the  magnetic  character  of  a  rock  can  be 
considered  to  be  made  up  of  residual  and  induced  magnetiza¬ 
tion.  Induced  magnetization  is  that  which  is  acquired  by  a 
rock  due  to  induction  by  the  earth's  present  magnetic  field. 
Residual  or  permanent  magnetization  is  that  acquired  by  the 
rock  as  it  cooled,  during  its  molten  state,  through  its 
Curie  temperature.  Such  magnetization  is  termed  permanent 
because  it  is  highly  resistive  to  removal  or  reversal.  It 
is  frequently  appreciable  in  magnitude  and  may  not  be  in  the 
same  direction  as  the  induced  magnetization  associated  with 
the  rock. 

Thus  the  magnetic  character  of  the  crust  is 
essentially  superimposed  on  that  of  the  geomagnetic  field 


susceptibility  -  a  measurement  of  the  ability  of  a  substance 

to  be  magnetized  under  the  influence  of  a 
magnetic  field. 
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and  any  attempt  to  analyze  the  field  associated  with  a 
portion  of  the  crust  must  involve  the  removal  of  the  regional 
gradient  due  to  the  terrestrial  field. 

The  magnetic  character  of  the  crust  varies  later¬ 
ally  in  response  to  changes  in  content  of  ferromagnetic 
materials.  However,  in  relation  to  depth  below  the  surface 
two  other  factors,  temperature  and  to  a  certain  extent 
pressure,  govern  the  amount  of  magnetization.  The  tempera¬ 
ture  gradient  in  the  crust  is  approximately  30°  C/km,  so 
that  at  a  depth  of  about  25  km.,  the  Curie  point  for  iron, 
namely  750°  C.,  is  reached.  At  the  Curie  point  ferromagnetic 
materials  become  paramagnetic  and  essentially  lose  most  of 
their  detectable  magnetic  character.  The  depth  extent,  of 
magnetization  within  the  crust  is  thus  temperature  controlled; 
the  influence  of  pressure  in  raising  the  Curie  point  of  a 
material  has  not  yet  been  fully  determined  but  more  will  be 
said  about  this  in  the  next  section.. 

1 . 4  Review  of  Publications  on  Related  Topics 

In  1941,  Vacquier  and  Affleck  made  a  statistical 
study  of  some  126  local  magnetic  anomalies  in  an  effort  to 
estimate  the  average  depth  to  the  Curie  point  geotherm. 

They  assumed  the  crust  to  be  divided  into  three  distinct 
layers:  the  first  being  a  sedimentary  layer  of  from  zero 

to  ten  miles  in  thickness,  a  magnetic  granitic  layer  and 
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an  underlying  layer  of  non-magnetic  or  homogeneously  magnetic 
material.  Using  the  narrow  anomalies  to  calculate  the  depth 
to  the  top  of  the  granitic  layer,  curves,  for  intensity  of 
magnetization  versus  this  depth  with  depth  to  the  bottom  of 
granitic  layer  or  to  the  bottom  of  the  layer  of  inhomogeneous 
magnetic  material  as  a  parameter,  were  constructed.  The 
data  were  then  fitted  to  these  curves  with  the  best  fit  result 
giving  an  average  depth  to  the  base  of  the  inhomogeneities 
of  12.5  miles.  Now  density  contrasts  at  much  greater  depths 
are  indicated  by  gravity  data  and  therefore  inhomogeneity  in 
magnetic  content  should  occur  if  the  material  at  these  depths 
can  be  magnetized.  However,  Vacquier  and  Affleck  found  no 
anomalies  which  indicated  a  source  at  a  greater  depth  and 
therefore  they  concluded  that  the  base  of  the  inhomogeneous 
layer  coincided  with  the  Curie  point  geotherm. 

The  above  estimate  of  the  depth  falls  within  the 
range  of  those  values  predicted  for  the  depth  to  the  bottom 
of  the  granitic  layer  by  Benfield  (19^0) .  His  work  was 
based  on  a  mathematical  treatment  of  heat  flow  measurements. 
However,  the  base  of  the  inhomogeneous  magnetic  layer  which 
Vacquier  and  Affleck  consider  to  be  the  Curie  point  geotherm 
could  be  the  Conrad  discontinuity.  The  depth  estimate  of 
12.5  miles  is  within  the  range  of  the  depth  to  the  Conrad  of 
5-20  miles  (Gutenberg,  1955) . 

In  1955,  Birch  published  a  paper  in  which  he 
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commented  on  the  work  of  Vacquler  and  Affleck  mentioned  above. 
He  pointed  out  that  they  may  have  excluded  broader  anomalies 
resulting  from  sources  at  greater  depths.  His  scepticism 
with  respect  to  their  results  seems  to  result  from  his  esti¬ 
mates  of  temperature  versus  depth  within  an  average  crust. 

For  various  models  of  crustal  structure  embodying  varied 
radioactive  content  and  heat  flow  contribution  from  below  the 
crust  and  also  due  to  radioactive  decay  Birch  constructed 
temperature  -  depth  curves  which  approach  the  Curie  point 
for  magnetite  (578°C)  in  the  30  to  40  km.  range.,  which  is 
from  1  1/2  to  2  times  greater  than  the  estimate  of  Vacquier 
and  Affleck.  Birch  also  states  that  ignorance  of  the  magnetic 
character  of  the  composition  of  the  rocks  at  these  depths 
raises  the  question  as  to  whether  further  analysis  of  regional 
magnetic  anomalies  might  not  shed  some  light  on  this  subject. 

Serson  and  Hannaford  (1957)  published  a  paper  in 
which  they  gave  their  results  from  a  statistical  analysis  of 
aeromagnetic  profiles  flown  over  continental  and  oceanic 
regions  in  1955*  They  calculated  autocovariance  functions 
for  profiles  flown  over  Western  Canada  and  the  Atlantic 
Ocean  and  then  fitted  models  to  these  curves.  They  concluded 
that  oceanic  magnetic  character  was  made  up  of  poles  at  a 
depth  of  one  kilometre  below  the  ocean  floor  with  an  intensity 
(susceptibility)  of  at  least  0.005  c.g.s.  units.  Furthermore, 
if  these  poles  had  corresponding  poles  their  depth  has  to  be 
greater  than  30  km.  (that  is  a  layer  of  dipoles  ^ 


30  km. 
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In  thickness) .  Under  the  continents,  they  said  the  rocks 
must  have  a  susceptibility  of  0.05  c.g.s.  or  greater  (for 
Western  Canada)  and  that  the  magnetic  structure  consisted 
of  a  thin  layer  at  a  depth  of  10  to  12  km.  or  a  layer  from 
4  1cm.  downwards. 

Perhaps  the  most  enlightening  contribution  on 
the  use  of  aeromagnetic  profiles  to  probe  into  the  earth’s 
magnetic  character  was  published  by  Alldredge  and  Van  Voorhis 
in  1961.  Using  profiles  with  a  length  of  at  least  2000 
miles  and  no  harmful  heading  changes  they  removed  that  part 
of  the  field  due  to  an  axially  centred  dipole  and  then 
plotted  a  histogram  of  crossover  distances.  This,  in  effect, 
was  a  very  simple  type  of  harmonic  analysis.  The  histogram 
showed  two  definite  ranges  for  crossover  distances,  those 
with  a  range  of  0  to  260  nautical  miles  and  those  in  the 
range  2100  to  5200  nautical  miles.  They  attributed  the 
former  to  crustal  magnetic  sources  and  the  latter  to  sources 
near  the  core-mantle  boundary.  The  large  gap  between  the 
crossover  distance  ranges  seems  to  indicate  that  the  mantle 
is  non-magnetic  in  character.  As  pointed  out  by  Serson  and 
Hannaford  also.  Alldredge  and  Van  Vorrhis  found  that  the 
majority  of  anomalies  occured  in  the  short-wavelength 
category.  In  their  case  they  found  that  93  per  cent  had  a 
crossover  distance  of  60  nautical  miles  or  less. 
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1 . 5  Object  of  the  Thesis 

The  purpose  of  this  thesis,  in  relation  to  the 
statement  of  the  problem  as  outlined  in  section  1.1,  is 
three-fold : 

(1)  To  gain  an  insight  intocontrasting  or  similar 
features  of  the  magnetic  character  of  continental  and 
oceanic  crusts  by  means  of  power  spectral  analysis  on 
aeromagnetic  data. 

(2)  To  show  that  aumerical  filtering  can  be  effectively 
applied  to  geophysical  data  to  reveal  information  which 
was  hitherto  inaccessible  through  ordinary  interpretational 
procedures . 

(3)  To  show  that  certain  anomalies  revealed  by  numerical 
filtering  can  be  interpreted  as  anomalous  Curie  point  geotherm 
structure  and  thus  perhaps  outline  a  practical  method  for 
detailed  analysis  of  magnetic  features  due  to  deep  sources. 
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II. 


POWER  SPECTRA 


2 . 1  Introduction 

Alldredge  and  Van  Voorhis  (1961)  found  that  when  they 
performed  a  harmonic  analysis  on  global  magnetic  data  the 
resulting  spectrum  revealed  two  possible  sources  of  terrestrial 
magnetic  anomalies.  One  was  the  core-mantle  boundary  with 
which  the  dipole  field  was  associated  and  the  other  was  the 
magnetic  material  in  the  crust. 

If  the  dipole  field,  which  is  composed  of  long 
period  signals,  were  removed,  then  the  shorter  period  crustal 
anomalies  could  be  studied  more  reliably  with  respect  to  a 
base  line  from  which  the  long  period  trend  had  been  removed. 

It  has  been  speculated  that  the  anomalies  associated  with  the 
Curie  point  geotherm,  if  they  exist,  would  be  broad  and  com¬ 
paratively  small  in  amplitude  due  to  the  depth  of  the  geotherm 
(approximately  10-60  miles  below  the  surface;  Serson  and 
Hannaford,  1957).  Because  of  these  two  characteristics,  the 
anomalies  are  likely  to  be  slightly  masked  by  those  of  near¬ 
surface  magnetic  materials.  Nevertheless,  the  energy  associ¬ 
ated  with  these  broad  anomalies  ought  to  be  comparatively 
greater  than  that  of  the  narrower  near  -  surface  structures. 

For  this  reason  and  because  of  the  fact  that  a  power  spectrum 
is  also  a  function  of  the  frequency  of  occurrence  of  a  certain 
period,  a  power  spectrum  analysis  was  performed  on  the  magnetic 
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data  available.  It  was  hoped  that  this  analysis  might  provide 
a  clearer  picture  of  the  type  of  anomalies  present  In  a 
typical  magnetic  profile,  as  well  as  perhaps,  to  define  more 
clearly  some  long  period  anomalies  whose  origins  might  be 
attributed  to  the  Curie  point  geotherm. 

Because  the  mathematical  approach  to  spectrum 
analysis  Is  generally  based  on  communications  engineering, 
one  seems  to  associate  such  an  analysis  with  a  stationary 
time  series,  that  Is  one  in  which  frequencies  are  expressed 
In  cycles/unit  time.  However,  the  spectrum  analysis  is 
equally  applicable  to  a  stationary  space  series  and  one  need 
only  substitute  the  space  variable  for  the  time  variable, 
thereby,  for  example,  obtaining  frequencies  in  cycles/unit 
length . 

If  a  measured  signal1  containing  noise  can  be 
described  as  an  approximately  stationary  random  process 
then  the  associated  autocovariance  function  or  power  spec¬ 
trum  will  generally  be  the  most  useful  display  for  determining 
the  presence  of  periodicities,  according  to  Blackman  and 
Tukey  (1958). 

The  power  spectrum  is  nothing  more  than  the 
Fourier  transform  of  the  autocovariance  function;  however. 


signal1  -  The  signal,  as  referred  to  in  this  thesis,  shall  be 
defined  as  a  message  which  is  recorded  on  a  strip 
chart  and  upon  which  is  superimposed  noise,  both 
instrumental  and  near-surface  geological  noise. 
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when  these  functions  must  be  determined  by  measurement  and 
computation  and  then  interpreted,  the  power  spectrum  lends 
itself  to  application  more  readily  than  the  autocovariance 
function.  The  reason  is  that  the  criteria  developed  so  far 
for  the  application  of  autocovariance  analysis  are  so  com¬ 
plicated  that  it  is  very  difficult  to  apply  them  in  practice. 

2.2  Definitions 

For  an  ideal  stationary  Gaussian  process  the  auto¬ 
covariance  function  is  defined  by 

X  /2 

C(X)  =  lim  ^  /  H(x)  •  H(x+X)  dx 

X-*°o  A  -X/2 

where  x  is  the  space  variable 

and  X  is  the  lag  in  units  of  length. 

C(X)  may  be  rewritten  as 

C (X)  =  f  °°  P(f)  •  ei27rfX  df 

-CO 

where  P(f)  =  lim  -  i  r  X/2  H(x)  •  e‘i27rfx  dx  i  2 
X-  X  '  J-X/2  ' 

and  P(f)  is  the  power  spectrum  of  the  process  considered. 

The  power  spectrum  can  also  be  expressed  as  the 
Fourier  transform  of  the  autocovariance  function,  that  is 

p(f)  =  C ( X )  e'i27rfX  dX 
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Since  both  P(f)  and  C(A)  are  even  functions  and 
because  only  positive  frequencies  are  of  Interest  we  have 
that 

P(f)  =  2/o°°  C(  A)  cos  (27rfA)  dA 

and  C ( A )  =  2 /“  P(f)  cos  (2irfA)  df 

2 . 3  Derivation  of  Power  Spectra  for  Equally  Spaced  Records 

In  practice  we  can  obtain  only  a  limited  number  of 
pieces  of  H(x)  of  finite  length;  each  piece  can  be  regarded 
as  a  sample  of  the  population  of  pieces  H(x)  of  the  same 
length.  The  reduction  of  such  data  will  therefore  only 
produce  estimates  of  the  power  density  spectrum  and  these 
will  be  subject  to  sampling  variations  and  to  biases  in  the 
usual  statistical  sense.  The  latter  statement  is  especially 
true  of  records  for  which  H(x)  is  obtainable  only  at  uniformly 
spaced  values  of  x  within  each  piece  H(x) .  Theoretical 
study  of  sampling  variability  and  bias  shows,  however,  that 
estimates  of  power  within  two  narrow  frequency  bands  separated 
by  a  small,  fixed  separation  will  become  more  precise  and  will 
approach  statistical  independence  with  respect  to  each  other 
as  we  consider  longer  and  longer  records. 

One  of  the  inherent  difficulties  in  using  equally 
spaced  data  is  that  of  aliasing  of  frequencies.  Suppose 
H(x)  is  available  only  for  equally  spaced  values  of  x,  say 
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x  =  0,  Ax,  2 Ax,  . .  nAx . 

then  it  is  clear  that  the  highest  frequency  which  can  be 
defined  by  such  a  sampling  interval  is 


where  f^  is  the  folding  (or  Nyquist)  frequency.  Then  the 
frequencies  f,  2f^  -  f,  2f^  +  f,  4f^  -  f,  4f^  +  f,  and  so  on, 
are  aliases  of  one  another  where  f  is  the  principal  alias. 

However,  with  the  data  used  for  this  thesis  the  problem  of 
aliasing  can  be  neglected  because  of  the  altitude  at  which 
measurements  were  made.  The  sampling  interval  employed  was 
three  miles  and  the  folding  frequency  is  therefore  l/6  cycles/mile 
which  is  undoubtedly  larger  than  the  highest  detectable  fre¬ 
quency  at  the  altitudes  at  which  the  records  were  made.'1'  The 
fortunate  exclusion  of  the  problem  of  aliasing  makes  the 
spectral  estimates  more  reliable  than  they  would  be  ordinarily. 

We  proceed  now  to  outline  a  method  for  power  spec¬ 
trum  estimation  from  a  uniformly  spaced  discrete  space  series 
of  finite  length;  this  outline  is  derived  from  the  more  com¬ 
plete  discussion  given  by  Blackman  and  Tukey. 

First,  let  H^,  H^, . .  Hn  be  the  space  series 

with  a  sampling  interval  Ax. 


1 


To  see  this,  the  reader  is  referred  to  the  segment  of 
actual  record  shown  in  figure  4. 
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Then  the  mean  lagged  products,  for  lag  Interval  AA  =  pAx, 
are  given  by 


1  q=n-pr 

n  -  pr  q=0  Hq  Hq+pr 


for  r  =  0,  1,  .  ,.,m  with  m  ^  n/p.  (In  all  cases  In  this 

thesis  p  was  equal  to  unity) .  Second,  raw  spectral  density 
estimates  are  computed  according  to  the  formula 


q=m-l 

S  =  AA  [  C  +  2  2  C  cos  +  c  COs  rrr] 

r  o  q=l  q  mm 


Third,  refined  spectral  density  estimates  are  computed  accord¬ 
ing  to 


U. 


=  a 
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00 

j=i 


ij 


[S  ^  .  + 
r+j 


S  .] 
r-  j 


The  a.  .  *s  are  the  coefficients  of  the  predetermined  spectral 
J 

window.  The  criteria  for  the  selection  of  an  appropriate 

spectral  window  for  a  particular  analysis  are  outlined  in 

Blackman  and  Tukey;  in  the  analysis  used  on  the  data  for  this 

thesis  the  coeffieients  for  the  spectral  window  chosen  were 

a 2q  =  0.50  and  a^  =  0.25  (all  others  being  zero) 

hence  we  may  write  U  =  0.25  S  ,  +  0.50  S  +0.25S 

0  r  r-1  r  r+1 

and  since  only  positive  frequencies  are  considered  the 

above  estimates  are  doubled. 


1 


A  discussion  of  the  spectral  window  employed  in  the  analysis 
and  its  function  can  be  found  in  appendix  1. 
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The  above  procedure  is  readily  adaptable  for  a 
high-speed  digital  computer  so  that  power  spectra  can  be 
obtained  with  a  minimum  of  computation  and  time  on  the  part 
of  the  interpreter. 

2 . 4  Power  Spectra  for  Digitized  Magnetic  Data 
A.  The  Data 

The  data  used  in  this  analysis  were  taken  from 
several  continuous  strip  chart  recordings  of  the  vertical 
component  of  the  magnetic  field  as  measured  during  the  air¬ 
borne  magnetic  surveys  (1953-1960)  of  the  Dominion  Observatory. 

The  instrument  used  in  making  the  measurements  was 
a  three-component  airborne  magnetometer  designed  and 
built  at  the  Dominion  Observatory  ( Serson.,  Mack,  Whitham; 

1957).  It  consisted  of  a  field-measuring  device  mounted  on 
a  horizontal  gyro-stabalized  platform  inside  the  aircraft. 

The  accuracy  with  which  the  instrument  measured  the  magnetic 
field  at  the  magnetometer  head  was  0.1  degrees  in  declination 
and  10  gammas  in  the  horizontal  and  vertical  components. 
However.,  its  accuracy  in  measuring  the  geomagnetic  field  was 
considerably  lower  because  of  uncertainties  in  the  corrections 
applied  to  the  readings  for  the  magnetic  field  of  the  air¬ 
craft  . 

The  error  associated  with  one  particular  observa¬ 
tion  made  in  flight  by  an  instrument  of  the  type  described 
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may  be  considered  to  be  the  sum  of  six  independent  errors 
which  are  listed  below: 

(a)  errors  in  measuring  the  magnetic  field  of  the 
magnetometer  with  respect  to  the  direction  reference  system, 

(b)  errors  in  the  direction  reference  system, 

(c)  errors  due  to  changes  in  the  magnetic  field  of 


the  aircraft, 

(d)  errors  due  to  magnetic  disturbances, 

(e)  errors  in  geographical  position, 

(f)  noise  contributed  by  the  measuring  equipment. 

Experimental  results,  obtained  by  flying  on 

different  headings  over  a  region  where  the  gradients  of  the 
magnetic  field  are  small  and  accurate  corrections  for  magnetic 
disturbances  are  available,  showed  that  the  probable  error 
of  a  single  observation,  after  correction  for  the  magnetic 
field  of  the  aircraft,  was  40  gammas  in  the  intensity  of  the 
horizontal  magnetic  vector  and  10  gammas  in  the  vertical 
component.  The  magnetic  field  of  the  aircraft  changed  over 
periods  of  a  few  weeks,  producing  an  uncertainty  in  the 
corrections  for  the  aircraft  field  which  increased  the  prob¬ 
able  error  of  survey  observations  to  60  gammas  in  the  hori¬ 
zontal  vector  and  30  gammas  in  the  vertical  component. 

Errors  in  navigation  and  plotting  the  flight-lines  and  the 
effect  of  magnetic  disturbances  resulted  in  a  further  increase 
in  the  probable  error  of  an  observation,  as  plotted  on  the 
charts,  to  100  gammas. 
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It  was  found  by  analysis  of  observations  made  in 
1955  that  the  chief  source  of  error  was  that  of  uncertainty 
in  geographical  position  of  the  order  of  4  miles  which  included 
errors  in  plotting  the  flight-lines  as  well  as  actual  errors 
in  navigation.  The  reason  for  such  inaccuracy  was  the  method 
of  navigation  which  was  by  map  reading,  with  4  or  5  location 
estimates  being  made  per  hour.  Aerial  photography  would  have 
increased  the  accuracy  but  also  increased  the  time  required 
for  plotting  the  results  of  the  survey.  Since  the  survey 
was  conducted  for  the  production  of  magnetic  charts  with  a 
scale  of  100  miles  to  the  inch,  the  accuracy  was  considered 
sufficient.  It  should  be  pointed  out  also  that  over  much  of 
Canada  radio  aids  to  navigation  such  as  Loran  or  Decca  are 
not  available;  however,  on  oceanic  flights  the  radio  aid 
Shoran  could  be  used  for  navigation. 

The  magnetometer  head  contained  three  orthogonal 
magnetic  detectors  of  the  saturated  transformer  type,  which 
produced  direct  currents  proportional  to  the  fore-and-aft 
and  transverse  horizontal  components  and  the  vertical  compon¬ 
ent.  These  and  the  heading  of  the  aircraft  were  then  fed 
into  an  analog  computer  which  displayed  continuously  on  strip 
chart  recordings  the  declination  in  degrees,  and  the  horizon¬ 
tal  and  vertical  field  components  in  gauss.  An  alternative 
display  automatically  presented  the  average  values  of  these 
quantities  over  successive  five  minute  intervals. 
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The  chart  speed  was  3/^  inch/minute  and  the  average 
value  of  the  field  was  taken  every  five  minutes  at  02  1/2, 

07  1/2,  12  l/2, .  minutes.  The  centre  of  the  averaging 

period  was  therefore  at  00,  05>  10,  ....  minutes  and  the 

average  value  of  the  field  during  each  five-minute  interval 
was  used  as  a  base  line  for  the  next  interval.  The  meter 
amplifier  was  zeroed  at  the  end  of  each  averaging  period  for 
about  ten  seconds  and  the  location  of  the  base  line  on  the 
chart  was  found  by  observing  the  position  of  the  trace  during 
this  time.  This  base  line  was  always  within  t  10  gammas  of 
the  zero  of  the  chart. 

The  chart  record  gives  the  value  of  the  field  at  the 
detector  head  and  this,  corrected  for  aircraft  gives  the 
estimated  true  value  of  the  field  measured  at  the  aircraft's 
altitude.  When  the  aircraft's  magnetic  heading  and  altitude 
are  fairly  constant  the  corrections  to  be  applied  remain  the 
same.  Care  was  taken  in  this  work  to  avoid  using  those  charts 
in  which  a  large  heading  change  was  detected  and  with  one 
exception  the  altitudes  at  which  the  field  was  recorded  were 
all  7000  feet  or  greater. 

The  aircraft  velocity  was  maintained  at  180  miles/ 
hour  as  closely  as  possible  which  meant  that  a  5  minute  aver¬ 
age  corresponded  to  15  miles  of  flight  path  and  that  one 
minute  or  3/^-  inch  of  chart  record  was  equivalent  to  3  miles 
of  flight  path.  Because  of  requirements  due  to  digitizing. 
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the  sampling  interval  was  chosen  to  be  3  miles.  The  location 
of  any  particular  anomaly  in  the  oceanic  profiles  will  no 
doubt  be  inaccurate  by  several  miles;  however,  the  features 
to  be  studied  are  of  such  a  scale  that  their  exact  locations 
are  not  a  binding  factor. 

Four  charts  for  the  continental  region  were  chosen, 
one  from  the  northern  Canadian  Shield  area  near  the  Arctic 
Circle  and  three  more  from  the  region  between  Ottawa  and 
Edmonton.  Four  charts  were  also  chosen  for  the  oceanic  regions, 
two  from  the  Pacific  and  two  from  the  Atlantic.  The  locations 
of  the  flight-lines  for  these  charts  are  shown  in  figures  1, 

2,  and  3. 

These  recordings  were  then  digitized  with  the  samp¬ 
ling  interval  given  above  and  the  resulting  data  were  plotted 
for  each  chart.  In  order  to  remove  the  trend  in  each  profile 
which  was  largely  due  to  the  dipole  field  of  the  earth  a 
base  line^  was  drawn  along  the  apparent  trend  through  the 
data  which  had  been  joined  together  by  a  continuous  curve. 

It  can  be  seen  from  the  resulting  power  spectra  that  this 
method  satisfactorily  removed  the  longer  periods  in  the  data 
due  to  the  dipole  field  which  was  of  no  interest  in  this 
analysis.  The  plotted  points  were  then  revalued  with  respect 
to  this  new  base  line  and  the  resulting  data  were  then  punched 
on  cards  and  made  ready  for  the  power  spectral  analysis  pro¬ 
gram  for  the  digital  computer. 

1  The  base  line  shape  was  almost  linear  in  the  case  of  oceanic 

profiles,  but  It  was  more  curvilinear  for  the  contin¬ 
ental  profiles.  See  figure  4  for  illustration. 
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B.  Reliability  of  Power  Spectral  Estimates 

The  reliability  of  the  spectral  estimates  was 
investigated  under  the  assumption  that  the  true  magnetic 
message  was  superimposed  on  Gaussian  noise,,  that  is,  noise 
which  exhibits  a  Gaussian  distribution. 

If  y^,  y^, . y^  are  independently  distributed 

according  to  a  Gaussian  distribution  with  a  zero  average  and 
a  unit  variance  then  the  quantity 

.  2  2,2,  ,2 

=  y1  +  y2  +  • • • •  +  yk 

follows  a  chi-square  distribution  with  k  degrees  of  freedom. 

2 

Now  defining  the  coefficient  of  variation  of  x^.  to  be  the 

2 

ratio  of  the  RMS  deviation  to  the  average  value  of  x^  we 

obtain  for  its  value  (2/k)1^  where  k  is  the  number  of  degrees 

of  freedom,  RMS  deviation  (X-^)  =  •\Aar  (x^)  =V  2k,  and 

average  (x^.  )  =  k.  Then  as  k  increases,  x^  becomes  less 

variable  and  furthermore,  the  above  relation  holds  for  any 
2 

multiple  of  x^  • 

Now  a  convenient  description  for  the  stability  or 

reliability  of  any  positive  estimate  is  its  equivalent  number 

of  degrees  of  freedom.  Since  the  above  relation  holds  for 

2 

any  multiple  of  x^  ^  confidence  intervals  may  be  set  up  for 

2 

the  estimates  based  on  the  above  consideration  of  the  xk 
distribution .  For  convenience  in  their  manipulation  in 
connection  with  power  spectra  these  intervals  are  usually 
expressed  in  decibels. 
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For  example,  if  an  observed  estimate  of  the  true 
average  power  associated  with  a  particular  frequency  band  is 
given  by  z  decibels,  then  provided  it  has  l8  equivalent 
degrees  of  freedom,  we  can  be  90%)  confident  that  the  true 
average  power  lies  within  a  5  db  interval  symmetrical  about 
the  value  z.  As  would  be  expected  the  smaller  the  decibel 
spread  the  greater  the  number  of  degrees  of  freedom  required 
to  assure  that  the  true  value  lies  within  that  particular 
interval.  Thus  the  reliability  of  a  particular  observed 
estimate  is  functionally  related  to  the  number  of  equivalent 
degrees  of  freedom  associated  with  that  estimate.  For  a  more 
detailed  description  of  these  decibel  spreads  and  their  deri¬ 
vation  see  Blackman  and  Tukey,  pages  21  -  23. 

The  number  of  degrees  of  freedom  will  now  be  shown 
to  be  a  function  of  the  amount  of  data  available  such  that  the 
greater  the  length  of  record  the  more  reliable  the  spectral 
estimates  will  be.  It  can  be  shown  that  for  practical  pur¬ 
poses  the  effective  total  length  of  record  may  be  expressed 
by  the  approximation 

~  Xn  -  Xm  (See  Blackman  and  Tukey  p .  17 ) 

where  =  effective  length  of  record 

Xn  =  total  length  of  record 

p  =  number  of  pieces  of  record. 

Xm  =  maximum  lag  distance  used. 

The  number  of  degrees  of  freedom  is  also  expressible  as  the 

number  of  elementary  frequency  bands  which,  for  design  pur- 


- 
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poses ,  may  be  approximated  by 


TT~  (See  Blackman  and  Tukey,  p.  24) 
Am 


In  order  to  estimate  the  amount  of  data  required 


for  a  known  number  of  degrees  of  freedom  one  must  consider 
the  degree  of  resolution  required  in  the  power  spectrum.  As 
stated  previously,  the  smoothing  process  for  spectral  estimates 
or  the  spectral  window  employed  was  that  with  coefficients 
0.25,  0.50,  0.25,  with  all  others  being  zero.  This  process 
supplies  estimates  every  1/2X  cycles/mile  where  Xm  is  defined 
to  be  the  maximum  lag  distance.  Now  adjacent  spectral  esti¬ 
mates  will  have  windows  which  overlap  each  other  considerably 
and  hence  they  will  not  be  completely  resolved.  However, 
consider  the  next- to-ad jacent  estimates,  the  overlap  will  be 
small  and  resolution  can  be  said  to  be  complete.  Hence 
resolution  may  be  defined  as 


Resolution  (in  cycles/mile) 


1 


In  addition,  as  shown  before,  the  Nyquist  or  folding 


frequency  is  given  as 


f 


-  where  Ax  is  the  sampling  interval. 

2Ax 


Blackman  and  Tukey  recommend  that  f^  be  chosen  such  that 

f„  =  3/2  f 

N  -v  max 


max 


where  f  is  the  maximum  frequency  one  wishes  to  analyze  for 
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in  the  spectrum  analysis. 

It  follows,  therefore,  that 


Ax  = 


1 


3f 


max 


and  the  necessary  number  of  lags  are  therefore  given  by 


m  = 


X 


m 


Ax 


3  X  f 
m  max 


The  necessary  number  of  data  points  then  follow  from 


n  = 


Now  Xm  f 

m  max 


in  =  (x '  +  I  x  )  3f 

v  n  3  nr  max 


=  3  X  1  f  +  p  X  f 

n  max  ^  m  max 

=  (3  ^  +  p)  (X  f  ) 

v  J  Xm  m  max' 


■^max 


maximum  frequency 
resolution 


=  number  of  resolved  bands. 


By  substitution, 

n  =  (1.5  k  +  p)  (number  of  resolved  bands). 

In  the  spectral  analysis  carried  out  on  the  aero- 
magnetic  data  it  was  found  that  a  resolution  of  1/300  cycles/ 
mile  was  required  in  order  to  define  the  peaks  in  the  spectrum 
clearly.  Furthermore,  it  was  ascertained  that  the  maximum 
frequency,  f  ,  required  would  be  of  the  order  of  1/10 
cycles  /mile  and  appropriately  the  folding  frequency  was 
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chosen  as  1/6  cycles/mile  which  is  approximately  3/2  f  . 

The  sampling  interval  defined  by  the  selection  of  the  folding 
frequency  is  given  by 

Ax  =  — —  =  3  miles 

2fN 


Then  the  number  of  lags  required  by  the  resolution  is 


_  Xm  3QQ 

TSx  3 


100 


The  number  of  resolved  bands  is  given  by 


X  f 
m  max 


=  300  x  1/9  = 


(since  f^  was  taken  as  1/6  cycles/mile,, 


f 


'max 


=  2/3  ‘  1/6  =  1/9  - 


33.3 

1/10  cycles/mile) 


As  pointed  out  previously.,  the  number  of  degrees  of 
freedom  required  for  a  decibel  spread  of  5  db  is  l8  and  there 
were  4  pieces  of  record  for  both  continental  and  oceanic 
regions.  Hence  the  total  number  of  data  points  required  on 
the  4  pieces  of  record  is  given  by 


n  =  (1.5  •  18  +  4)  33.3  ~  1000 

The  continental  region  spectrum  was  calculated 
using  ll6l  points  and  that  for  the  oceanic  region  using  1257 
points.  The  noise  level  on  each  spectrum  is  estimated  to  be 
about  50  units  of  average  power  or  70  db .  The  90/  confidence 
range  about  this  level  is  5  db  which  is  symmetrical  about 
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POWER  SPECTRA  FOR  CONTINENTAL  PROFILES  FIGURE 
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POWER  SPECTRA  FOR  OCEANIC  PROFILES  FIGURE 
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the  70  db  level.  That  Is,  there  Is  only  a  5 $  probability  that 
a  peak  in  the  spectrum,  which  is  located  above  the  72.5  db 
level  (53  power  units),  will  be  due  to  noise  in  the  magnetic 
signal.  Therefore  any  peak  found  above  the  53  unit  level  can 
be  taken  as  statistically  significant  with  a  95$  confidence. 

C.  Results  of  Power  Spectral  Analysis 

The  digitized  data  for  each  flight-line  was  pro¬ 
cessed  using  the  power  spectral  analysis  program1  for  the 
high-speed  computer.  A  print-out  was  obtained  for  each 
flight-line  which  gave  the  power  associated  with  each  esti¬ 
mate  at  a  particular  frequency.  Since  the  Nyquist  frequency 
was  taken  to  be  1/6  cycles/mile  then,  as  shown  before,  an 
estimate  was  obtained  every  1/600  cycles/mile  from  0  to  1/6 
cycles/mile . 

The  power  spectra  for  the  aeromagnetic  profiles 
flown  over  the  continental  and  oceanic  regions  are  shown  in 
figures  5  and  6.  The  spectral  estimates  associated  with  these 
spectra  are  not  statistically  significant  because  of  the 
insufficient  record  length  in  each  case.  However,  it  was 
assumed  that  each  set  of  four  profiles  represented  a  random 
sample  of  the  infinite  population  of  profiles  of  the  magnetic 
character  of  the  continental  and  oceanic  regions;  thus,  in 
order  to  obtain  a  statistically  significant  power  spectrum 


1This  program  was  written  for  the  IBM  1620  digital  computer 
by  R.  M.  Ellis,  Department  of  Physics,  University 
of  Alberta. 
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for  each  of  the  two  regions  the  spectral  estimates  for  each 
region  (continental  and  oceanic)  were  weighted  according  to 
the  number  of  data  points  in  the  particular  length  of  record. 
The  number  of  data  points  and  corresponding  weighting  factors 
for  each  of  the  flight-lines  involved  are  shown  in  table  1. 

Table  1 


Region 

Flight  No. 

No.  of  Data  Pts. 

Weighting  Factor 

Continental 

2 

180 

0.1550 

2 

378 

0.3256 

7 

369 

0.3178 

26 

234 

0.2016 

ll6l 

1.0000 

Oceanic 

10 

612 

0.4869 

19 

159 

0 . 1265 

20 

270 

0.2148 

23 

216 

0.1718 

1257  1.0000 

The  weighted  estimates  for  the  continental 
and  oceanic  regions  were  then  plotted  and  the  resulting 
spectra  are  shown  in  figures?  and  8. 
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WEIGHTED  CONTINENTAL  POWER  SPECTRUM  FIGURE 
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WEIGHTED  OCEANIC  POWER  SPECTRUM  FIGURE 
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In  interpreting  these  power  spectra  one  must  not 
assume  that,  because  the  spectra  contain  certain  dominant 
frequencies,  the  magnetic  character  of  the  crust  is  necessar¬ 
ily  periodic.  Such  a  conclusion  is  unwarranted  because  of 
observations  made  of  surface  geology;  however,  there  are 
certain  marked  similarities  and  dissimilarities  between  the 
spectra  for  the  two  regions  which  should  be  investigated.  A 
detailed  discussion  of  these  features  will  be  given  in  chapter 
4;  at  this  time  we  wish  only  to  outline  the  reasons  for  the 
application  of  numerical  filtering  to  the  data. 

As  mentioned  in  section  2.1,  it  was  a  general  belief 
that  anomalies  of  high  intensity,  due  to  near-surface  magnetic 
structures,  might  mask  certain  long  period  anomalies  with 
relatively  low  amplitudes.  In  order  to  reveal  these  long 
period  features  more  clearly  and  accurately  than  interpolation 
by  the  human  eye  could  accomplish,  it  was  decided  that  the 
data  should  be  convoluted  with  appropriate  discrete  filters 
or  weighting  functions.  This  type  of  smoothing  process  could 
also  be  applied  more  quickly  than  that  of  smoothing  by  hand. 

Although  the  peaks  on  both  spectra  appear  to  occur 
at  approximately  the  same  frequency  in  all  cases,  we  shall  in¬ 
vestigate  fully  only  one  significant,  common  peak,  that 
which  occurs  at  a  period  of  approximately  250  miles  in  both 
spectra.  Since  the  trend  due  to  the  dipole  field  was  removed. 
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this  long  period  is  evidently  due  to  some  broad,  deep  struc¬ 
ture  which  seems  to  be  common  to  both  continental  and  oceanic 
magnetic  crusts .  One  such  structure  which  might  cause  such 
broad  anomalies  could  be  an  upwarping  or  downwarping  of  the 
Curie  point  geotherm.  This  possibility  will  therefore  be 
investigated  in  the  next  chapter  with  the  aid  of  numerical 
filtering . 
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III. 


NUMERICAL  FILTERING 


3 . 1  Introduction 

The  procedure  of  applying  a  weighting  function  or 
filter  to  a  set  of  data  in  order  to  eliminate  undesirable 
noise  from  the  signal  is  nothing  more  than  a  type  of  smooth¬ 
ing  process.  Just  as  the  formation  of  the  autocovariance 
function  is  the  smoothing  of  the  signal  by  the  signal  itself, 
the  action  of  a  filter  (a  low-pass  or  bandpass  filter)  is  to 
remove  by  smoothing  undesirable  high  frequency  components. 

In  numerical  filtering  of  equi-spaced  data  the 
filter  form  is  that  of  a  weighting  function  which  is  composed 
of  several  discrete  points  separated  from  one  another  by  the 
sampling  interval  for  the  data  in  question.  The  shape  of  this 
weighting  function  is  governed  by  the  frequency  characteristics 
which  the  filter  must  have  in  order  to  perform  the  required 
operation . 

The  proper  application  of  a  filter  to  any  type  of 
data  inevitably  requires  some  predetermined  information 
about  the  message  to  be  extracted  from  the  data.  In  geophysics, 
such  a  priori  knowledge  might  be  supplied  through  the  geology 
or  as  in  this  case,  by  a  power  spectral  analysis,  provided 
the  input  signal  is  assumed  to  be  random.  Such  an  analysis 
defines  the  range  of  frequencies  against  which  the  filter  must 
discriminate . 


41 


,  a  a  6  c  ' ;  ' 

:  ■  •  fi  oc  -.V  eri :  ,  1  nr  :  . 

'■  -i.r  i  .  .  .  bo  l  g  B  uv ;  - 

■'  '  ■  '  a X  bb roGoorna  vci  bvobgb 

£  bI  p-L -  i bv/  s/rc-  cr&rfB  ar  jfub: 

■  '  »  ■  §  a;|  ...  1  ,  b  ctfrf;  .  s 

-  -i  :-v  •  "o  b  jd  ooiri'BolXqqB  rioqp{m  e nri‘ 

•  i  OCBB\SB  ,  BB  '  f  OJ  BBB' G  3  9IU  BBj  JCjBBG 

t  '  1  ■../'■■■  :  .  . 


..  :  • "  '  - ::  Gi:!.; 


42 


In  general,  there  are  three  basic  types  of  filters: 
low-pass,  bandpass,  and  high-pass.  The  low-pass  discriminates 
against  high  frequencies,  the  bandpass  against  both  low  and 
high  frequencies,  and  the  high-pass  against  low  frequencies. 

In  most  geological  problems  the  filter  required  will  either 
be  of  the  low-pass  or  bandpass  type  since  most  of  the  noise 
in  the  signal  is  usually  caused  by  near- surface  geologic  effects 
which  are  composed  mainly  of  high  frequency  components. 

One  distinct  advantage  of  numerical  over  electrical 
filters  is  that  the  weighting  function  to  be  applied  to  the 
data  does  not  need  to  be  physically  realizable.  That  is,  one 
does  not  have  to  reproduce  the  weighting  function  as  the  out¬ 
put  of  an  electrical  circuit.  Once  the  frequency  character¬ 
istics  of  the  weighting  function  have  been  determined  to  be 
those  required  to  perform  the  necessary  filtering  operation 
the  filter  is  thus  designed.  In  this  connection,  it  is  best 
to  point  out  that  linear  filters  behave  well  in  the  neigh¬ 
bourhood  of  the  signal  origin.  They  are,  in  fact,  independent 
of  signal  origin;  but,  it  is  best  to  choose  an  origin  which 
is  at  lease  one  weighting  function  away  from  the  approximate 
location  of  the  anomaly  being  sought  in  order  that  part  of 
the  desired  message  not  be  lost  in  filtering. 
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3 . 2  Filter  Design 

In  designing  a  filter  or  weighting  function  there 
are  several  criteria  to  which  one  must  adhere  in  order  that 
the  filter  perform  the  desired  operation  reliably.  These 
criteria  in  order  of  importance  are: 

(1)  the  filter  must  preserve  the  message  frequency  while 
discriminating  against  undesirable  high  frequency  noise. 

(2)  the  filter  must  not  distort  the  true  shape  and  ampli¬ 
tude  of  the  message. 

(3)  there  should  be  no  phase  shift  due  to  filter  application, 
if  at  all  possible. 

(4)  the  discrete  weighting  function  factors  should  be 
normalized  with  respect  to  unity  in  order  that  unit 
response  be  given  to  a  function  of  infinite  period. 

In  order  to  determine  if  a  filter  satisfies  the 
first  criterium  one  must  study  the  frequency  characteristics 
of  the  weighting  function.  This  is  usually  accomplished  by 
computing  the  Fourier  transform  of  the  weighting  function 
which,  in  effect,  gives  the  response  of  the  filter  to  a 
given  input  with  respect  to  the  frequency  of  that  input. 

In  this  connection,  it  is  best  to  point  out  that  one  can 
proceed  from  the  Fourier  transform  to  filter  shape  as 
Ormsby  describes  (i960),  that  is,  one  can  design  the  desired 
transform  function  and  then  derive  the  corresponding  weighting 
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function,  or  one  can  experiment  with  filter  shapes  directly 
(Goldstein,  1962) .  The  latter  method  is  quite  satisfactory 
for  low-pass  filters;  in  this  method  one  must  determine 
the  frequency  characteristics  of  the  filter  each  time  the 
weighting  function  shape  is  changed. 

It  was  found  that  this  procedure  was  adequate  if  one 
made  the  weighting  function  the  same  width  as  the  width,  (in 
miles)  of  the  message  sought.  This  period  was  taken  from  the 
power  spectrum  for  the  particular  flight  profile  to  which  the 
filter  was  to  be  applied.  In  most  cases  this  period  never 
varied  more  than  20  to  30  miles  on  either  side  of  the  average 
long  period  of  250  miles  appearing  on  both  continental  and 
oceanic  power  spectra.  The  shape  of  the  weighting  function 
was  then  altered  until  the  Power  Frequency  Transfer  Function"*- 
satisfied  the  specified  frequency  characteristics. 

The  definition  of  the  Transfer  Function  is  derived 
from  the  following  mathematical  analysis  of  the  filtering 
process . 

The  input  to  a  linear  filter  is  related  to  its 
output  by  the  following  integral  equation  known  as  the 
convolution  integral: 

^o  W  =  /  I  ljU-6)  W(x)  dx  - (1) 

where  x  =  space  variable 

W  =  weighting  function 
6  =  lag 


■'"Power  Frequency  Transfer  Function  -  the  absolute  value  of  the 

Fourier  transform  of  the  weighting  function,  all  squared. 


-  "i  :■  up- 

•  '!  W 

r":.-;  '  ■ op 

P  £  £ 

rl  ■'■:)  ■ .  y 

- 

J.B  !  :  ■■  .■• 


45 


To  express  a  function  of  a  time  or  space  variable 
as  a  function  of  frequency  one  need  only  take  the  Fourier 
transform  of  the  function.  For  example, 

E(a>)  =  -L  /  “  l(x)  e+lffix  dx  ....(2) 

»/27  J  -00 

and  inversion  gives, 

b  /  \  1  _  oo  -r-i  /  \  -iCLSX  J 

|(x)  =  — =  f  E(o))  e  da> 

v  2tt  -°° 

where  cd/2t r  is  the  frequency  in  cycles/mile. 

The  Fourier  transform  of  the  integral  equation  (l) 
gives  a  simple  input-output  relation  in  terms  of  the  frequency 
variable,  go. 


Eo(cjd)  =  E.(o))  Y(o))  . (3) 

where  the  functions  Eq(o)),  E^co),  Y(<u)  are  the  Fourier  trans¬ 
forms  of  the  functions  £,q(x),  |^(x),  W(x)  respectively. 

From  (3)  follows  the  definition  of  the  Power 
Frequency  Transfer  Function, 

Y(cd)  j  2  =  (-jh)2  . (4) 

The  analogy  between  the  above  expression  and  that  of  the  ratio 
of  output  over  input  voltage  in  an  electrical  circuit  is 
obvious.  Hence  the  name  Power  Frequency  Transfer  Function. 

“1  00 

Y(cd)  =  -  r  W(x)  [  cos  ( cox )  +  j  sin(cQx)]  dx  ...(5) 

/SF  J-“ 


Now 


-  ■  3-;-  p.  a  ,10  ao.t^oar/ij  &  as 

■  .3..,  •• :  -  f.  1)  t  .•.■■■■■: ;  lei-KJo  •  e:.  ' 

■  :■  WQiUrt  ic)  :  J'l-J. 

■ '  .oio  ...  /*?  'IP'., i .  PI  Yf-;OJ..p  v‘:i.a 

■  -  -  r  '  :  .  ' 

t  ujc  .iiYdo 

■n'p'v 


Since  the  weighting  function  is  given  by  a  series  of  discrete 
points  separated  from  one  another  by  the  sampling  interval 


Ax,  then. 


Y(cd)  =  —  2  W(x1)  [cos(o3Xji)  +  j  sin(cox1)  ]  ...(6) 

J2 TT  X±=-«» 


00 


For  a  symmetrical  weighting  function. 


Y(cq)  =  -  Ax  2  W(x.  )  [  cos  ( cox . )  +  j  sin  (cox  . )  ]  .  .  .  (7 ) 

7 r  Xj_=o  ii  i 


Let  Y(a>)  =  J2-  Ax  ( A  +  jB) 


.  .  |  Y( at)  |  =  Ax  |  7a2  +B2  |  . (8) 

V  7 T 

and  the  Power  Frequency  Transfer  Function  is  given  by, 

I  Y(cu)  |  2  =  21  — (A2  +  B2)  . (9) 

IT 

00 

where  A  =  2  W(x.)  cos  (cox.) 


00  /  \  /  \ 

B  =  .2  W(x .  )  sin(cux . ) 

Xj_=0  1  ! 

N 

and  2  W(x.)  =  1  for  a  weighting  function  of  length  NAx . 

i=l  1 


The  third  of  the  aforementioned  criteria  is  satisfied 


if  the  weighting  function  is  symmetrical  with  respect  to  its 
centre  point.  This  construction  also  simplifies  the  summation 
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given  in  equation  (7)  above. 

A  more  convenient  form  of  the  Transfer  Function 
for  plotting  is  obtained  if  one  converts  the  quantity  to 
decibel  units. 

From  equation  (9)  then, 

|  Y(<o)  |  2  =10  log...  |  Y(a>)  |  2  ...(10) 

db  1U 

It  should  be  mentioned  that  the  summations  in 
equation  (9)*  though  taken  for  from  0  to  00  really  only 
apply,  ef f ectively,  over  the  finite  length  of  the  weighting 
function  since  beyond  the  end  points  the  function  is  zero 
valued.  These  calculations  can  be  easily  programmed  for  the 
digital  computer  and  one  can  then  obtain  the  Transfer  Function 
quite  rapidly. 

One  of  the  weighting  functions  used  in  the  numer¬ 
ical  filtering  of  the  flight  data  is  shown  in  figure  9- 
This  filter  was  obtained  after  seven  refinements  of  the  orig¬ 
inal  trial  filter  and  satisfies  not  only  the  frequency  require¬ 
ments  but  also  criterium  2. 

The  Transfer  Function  for  this  filter  is  shown  in 
figure  10.  It  can  be  seen  that  the  filter  is  reasonably 
sharp  and  therefore  will  discriminate  against  high  frequency 
components  quite  adequately.  The  very  low  frequency  compon¬ 
ents  are  not  a  factor  here  since,  as  can  be  seen  from  the 
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TEST  FOR  DISTORTION  DUE  TO  FILTERING  FIGURE  11 
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power  spectra,  they  have  been  adequately  removed.  In  the 
frequency  range  of  Interest,  the  filter  is  very  slightly 
suppressive;  however.  It  Is  very  close  to  the  zero  decibel 
axis  and  therefore  will  have  approximately  a  1:1  ratio  for 
response/input.  This  factor  is  very  important  in  order  that 
criterium  2  be  satisfied  as  closely  as  possible. 

In  addition  to  this  factor,  the  weighting  function 
must  be  so  constructed  that  the  shape  of  the  message  is  not 
appreciably  distorted.  This  is  necessary  since  depth  calcula¬ 
tions  for  magnetic  anomalies  are  very  dependent  on  the  shape 
of  the  anomaly.  In  order  to  check  the  filter  out  on  this 
latter  stipulation  an  anomaly  was  constructed  with  test  data 
whose  length  was  identical  to  that  of  the  weighting  function. 
This  test  data  was  then  convoluted  with  the  filter  to  deter¬ 
mine  the  amount  of  distortion  due  to  filtering.  The  test 
data  anomaly  and  the  filtered  result  are  shown  in  figure  11. 

The  filter  thus  duplicates  the  anomaly  slopes 
quite  accurately  and  the  amplitude  suppression  is  only  25 
gammas  in  300  or  8  per  cent.  The  filter  can  also  be  suit¬ 
ably  altered  in  shape  and  length  to  handle  messages  with 
periods  in  the  200  to  300  mile  range. 

Having  succeeded  in  designing  an  appropriate  filter 
the  various  flight  profiles  were  then  convoluted1  with  this 
filter  (for  various  periods)  producing  a  filtered  profile. 

An  example  of  such  a  profile  is  given  in  figure  12 . 

^The  mathematical  procedure  for  this  convolution  which  was 

programmed  for  the  digital  computer  is  outlined  in 
appendix  2. 
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EXAMPLE  OF  A  FILTERED  PROFILE  FIGURE  12 
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3 . 3  Results  of  Numerical  Filtering 

As  mentioned  in  section  3*2,  the  filter  produces 
its  best  accuracy  when  applied  to  a  signal  containing  an 
anomaly  whose  width  is  identical  to  that  of  the  filter  itself. 
Consequently,  when  the  filtered  data  were  obtained  from  the 
digital  computer  they  were  then  scanned  for  anomalous  features 
whose  widths  were  within  ±  3  miles  of  that  of  the  filter 
applied.  It  was  found  that  an  anomaly,  whose  width  was  greater 
than  or  less  than  that  of  the  filter  applied  by  more  than  3 
miles,  was  not  shortened  or  lengthened  as  a  result  of  filter 
application.  That  is  the  width  of  the  anomaly  was  relatively 
constant  under  filter  application  as  long  as  the  filter  width 
was  within  10  -  20  miles  of  that  of  the  anomaly.  Therefore, 
when  an  anomaly  of  greater  or  lesser  width  than  that  of  the 
filter  applied  was  found  the  appropriate  filter  of  proper 
width  was  then  applied  insuring  a  minimum  of  distortion  of 
anomaly  shape . 

Because  of  the  relatively  large  periods  of  the 
anomalies  being  sought  it  was  difficult  to  find  more  than  one 
on  any  particular  flight  profile.  Only  on  flights  2  and  10 
were  two  anomalies  found  in  each  case;  flights  J,  20,  and  23 
produced  one  each.  The  anomalies  produced  by  filtering 
profiles  for  flights  17  (Atlantic  ocean)  and  26  (continental) 
were  comparatively  small  in  width  (having  a  period  much  less 
than  200  miles)  and  so  did  not  qualify  for  interpretation. 
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In  general,  it  appears,  from  the  seven  anomalies 
analysed,  that  the  oceanic  anomalies  are  relatively  smaller 
in  amplitude  than  those  of  the  continental  profiles.  It 
should  also  be  pointed  out  that  all  of  the  oceanic  anomalies 
were  positive  while  those  of  the  continental  region  were  all 
negative.  It  is  suspected  that  this  latter  feature  may  not 
be  significant  in  view  of  the  very  small  number  of  anomalies 
found  on  the  profiles.  A.  general  description  of  the  anomalies 
found  is  given  in  table  2. 


Region 


Table 


A  n  o  m  a  1  \ 


«y 


Features 


Flight  No.  No.  Anomaly  Max.  Amplitude 

_  _  width  (miles)  ( gammas ) _ 


Oceanic 

10 

A 

117 

44 

B 

102 

107 

20 

A 

108 

67 

23 

A 

144 

228 

Continental 

2 

A 

100 

188 

B 

121 

348 

7 

A 

126 

204 

It 

should  be 

noted  that 

anomaly 

width  as 

referred 

to  in  table 

2  is  that 

portion  of 

the  anomaly  which 

occurs 

between  the 

zero  amplitude  crossings  near 

each  end 

of  the 

anomaly. 
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A  more  complete  discussion  of  these  anomalies  and 
their  interpretation  will  be  given  in  chapter  4.  Diagrams 
of  the  anomalies  will  be  given  there  also. 
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IV. 


INTERPRETATION  OF  RESULTS 


4 . 1  Interpretation  of  Power  Spectra 

An  interpretation  of  the  contrasting  and  the  similar 
features  of  the  weighted  power  spectra  shown  in  figures  7  and 
8  in  chapter  2  will  now  be  given. 

A  power  spectrum  for  continuous  or  discrete  data 
displays  estimates  of  the  true  average  power  associated  with 
particular  frequency  band  widths  over  a  specified  range  of 
frequencies.  The  amplitude  of  any  estimate  is  a  function  of 
the  frequency  of  occurrence  of  the  corresponding  period  as 
well  as  the  energy  associated  with  this  period. 

Through  observation  of  the  eight  power  spectra 
obtained  (see  figures  5  and  6,  chapter  2),  the  maximum  noise 
level  was  approximated  as  50  db,  due  to  the  gradual  decline 
from  this  level  on  all  of  the  spectra  in  the  shorter  period 
range  (<  30  miles) .  Thus  the  four  peaks  which  appear  on 
each  of  the  weighted  spectra  can  be  said  to  be  statistically 
significant  in  view  of  the  analysis  of  reliability  of  spectral 
estimates  given  in  section  2.4.  That  is,  one  can  be  95  per¬ 
cent  confident  that  these  peaks  are  not  due  to  some  anomalous 
feature  of  the  noise  spectrum. 

Bearing  these  facts  in  mind,  let  us  analyze  the 
frequency  characteristics  of  the  continental  and  oceanic 
spectra  shown  in  figures  7  and  8.  In  both  cases  there  is  a 


56 


i'  V  .  ■' 

■'  " ir/; o  L o c  1 

■  'io  ilq-u  enl?  .  •  jjtoxv  ••  :‘t 

i"W  l-  ’  ; '  ■/  ",  i ■:  nv.  '  ; ::  :  ;■  \t 

0  ,,  1  v-  ?  ^ :  99#)  Oi'ilS; 

™  -  :  0;,,;  ,;;:i  be:?  Am^o'xeq.:  9,sv/  love  I 

v  '  i  I  [/'  "‘o  l ..  v  -  ,  :ri  fnov 

v  ' •  :'v.. ' ''  v T  ■  .  (  3  o.C 9f  ::  ■ )  e:::  • , 

.  33  ”  :  .. r . ; r,o  9~EX.i'  ' 

' 

i ..  j  :9':  -e /•  ;p£  r ? 


57 


sharp  rise  from  the  zero  frequency  value  to  the  peak  at  the 
frequency  of  approximately  1/250  cycles/mile.  The  large  gap 
between  these  two  frequencies  is  prominent  in  both  spectra  and 
shows  that  there  is  only  a  comparatively  small  amount  of 
energy  associated  with  those  periods  which  lie  between  infin¬ 
ity  and  250  miles.  Thus  it  appears  that  the  removal  of  the 
dipole  field  trend  causes  those  long  periods,  supposedly 
associated  with  the  core-mantle  boundary,  to  disappear.  This 
result  seems  to  agree  with  that  of  Alldredge  and  Van  Voorhis 
(1961).  Recall  from  chapter  1  that,  using  a  simple  type  of 
harmonic  analysis,  they  attributed  those  features  of  extremely 
long  period  to  the  core-mantle  boundary  sources  and  those 
whose  crossover  distances  lay  between  zero  and  260  nautical 
miles  were  said  to  be  due  to  crustal  magnetic  materials.  It 
appears  that  the  sources  which  produced  the  four  anomalous 
peaks  on  each  of  the  power  spectra  can  be  classed  as  crustal 
or  near-crustal  structures. 

Although  these  spectra  are  only  two  samples  from 
the  infinite  population  of  possible  spectra  for  crustal 
magnetic  character,  there  could  be  some  significance  to  the 
fact  that  the  location  of  each  of  the  four  peaks  on  the 
continental  spectrum  is  approximately  coincident  with  that 
of  the  corresponding  peak  on  the  oceanic  spectrum.  The  four 
periods  under  consideration  are  in  the  neighbourhood  of 
250,  65  *  50,  and  40  miles.  These  periods  correspond  to 
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anomaly  widths  of  roughtly  125 >  33 >  25,  and  20  miles  respec¬ 
tively,  since  an  anomaly  for  a  relatively  broad  structure 
has  a  shape  which  approximates  half  a  sine  wave. 

One  possible  explanation  Is  that  each  of  these 
periods  Is  characteristic  of  structures  which  occur  at  a 
certain  depth;  that  Is,  the  higher  frequency  anomalies  are 
due  to  near-surface  narrow  magnetic  structures  and  as  the 
depth  increases  the  width  of  the  source  must  increase  in  order 
to  compensate  for  attenuation  of  the  signal  sent  out  by  the 
source.  The  attenuation  of  the  signal  from  a  magnetic  source 
below  the  surface  is  proportional  to  both  depth  and  frequency 
(Dean,  1958) .  Hence  the  signal  from  a  narrow  structure 
(higher  frequency)  would  become  gradually  negligible  if  the 
depth  to  this  structure  were  increased.  Also,  in  most  areas, 
there  is  a  rapid  lateral  change  in  susceptibility  contrast 
among  the  various  rock  types  on  or  near  the  surface.  These 
changing  contrasts,  though  small,  can  give  rise  to  narrow 
anomalies  and  broad  anomalies  due  to  near-surface  structure 
are  not  detectable.  Due  to  attenuation  these  sharp  lateral 
changes  gradually  become  smoothed  over  and  in  order  to 
produce  an  anomaly  in  the  earth's  field  measured  at  the  sur¬ 
face  a  greater  susceptibility  contrast  is  required  as  the 
depth  increases.  This,  in  turn,  requires  a  greater  amount 
of  lateral  distance  along  the  structure  in  order  that  a 
sufficient  contrast  is  achieved,  thus  giving  rise  to  broader 
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anomalies.  Therefore,  one  might  say  that  the  period  of  a 
magnetic  anomaly  increases  with  the  depth  to  the  structure  in 
which  the  source  is  located. 

This  increase  in  anomaly  period  will  continue  until 
the  Curie  point  geotherm  is  reached.  At  this  depth  there 
will  be  a  sharp  drop  in  susceptibility  and  the  rocks  will 
essentially  cease  to  be  ferromagnetic  any  longer.  It  is 
therefore  hypothesized  that  this  sharp  drop  in  susceptibility 
corresponds  to  the  sharp  drop  from  the  250  mile  peak  to  the 
zero  frequency  value  on  each  of  the  power  spectra;  that  is, 
the  250  mile  period  anomalies  are  due  to  structures  associ¬ 
ated  with  the  Curie  point  geotherm. 

Having  analyzed  the  frequency  characteristics,  let 
us  examine  the  amplitudes  of  the  peaks  occurring  on  both 
spectra.  As  stated  previously,  in  more  general  terms,  the 
amount  of  power  associated  with  any  one  particular  spectral 
estimate  is  dependent  on  the  frequency  of  occurrence  of  the 
corresponding  period  and  on  the  amount  of  magnetic  intensity 
associated  with  the  structure  which  gives  rise  to  an  anomaly 
with  such  a  period.  It  is  evident  that  the  three  smaller 
peaks  on  both  spectra  occur  in  roughly  the  same  range  of 
amplitude  while  that  of  the  250  mile  period  is  much  larger 
for  the  continental  region  than  that  of  the  oceanic  region. 
The  frequency  of  occurrence  of  the  smaller  structures  associ¬ 
ated  with  the  smaller  periods  will  no  doubt  be  much  higher 
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than  that  of  the  structures  with  the  250  mile  period.  Because 
of  this  high  rate  of  occurrence  any  difference  in  magnetic 
intensity  could  be  counter-balanced  by  a  slight  difference  in 
frequency  of  occurrence  between  continental  and  oceanic  regions . 
However,  this  should  not  be  the  case  for  the  longer  period 
anomalies.  The  investigation  of  these  long  periods  by  numer¬ 
ical  filtering  showed  that  the  frequency  of  occurrence  for 
both  regions  was  roughly  equal  but  that  the  amplitudes  associ¬ 
ated  with  the  broad  anomalies  was  generally  less  for  the 
oceanic  region.  Although  these  amplitudes  are  not  completely 
accurate,  the  trend  is  still  evident  and  points  to  the  fact 
that  if  these  broad  anomalies  occur  at  approximately  the  same 
depth  in  both  continental  and  oceanic  crusts  then  the  intensity 
of  magnetization  is  greater  for  the  continental  structures 
than  for  the  oceanic  ones.  This  hypothesis  is  supported  by 
the  findings  of  Serson  and  Hannaford  (1957)- 

4 . 2  Interpretation  of  Numerically  Filtered  Anomalies 

In  section  3-3  a  brief  description  of  the  seven 
long  period  anomalies  revealed  by  numerical  filtering  was 
given.  In  this  section  a  more  complete  description  of  the 
anomalies  will  be  given  and  it  will  also  be  shown  that  the 
sources  of  these  anomalies  can  be  interpreted  as  anomalous 
Curie  point  geotherm  structures. 
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First  of  all,  one  must  realize  that,  although  the 
long  period  trends  (250  miles)  are  quite  prominent  on  the 
power  spectra  for  both  regions,  numerical  filtering  cannot 
reveal  these  broad  anomalies  with  the  same  clarity  as  a 
power  spectrum  might  detect  their  presence.  This  is  due  to 
the  fact  that  the  superposition  of  higher  frequency  near¬ 
surface  features  on  the  broad  anomalies  can  distort  them 
beyond  recognition  by  a  filtering  technique.  Hence  only  in 
those  areas  where  the  magnetic  character  of  shallower  structures 
is  relatively  homogeneous  can  the  broad  anomalies  be  extracted 
in  their  entirety  by  filtering.  This, coupled  with  the  fact 
that  the  anomaly  width  was  large  compared  to  the  flight 
lengths  available,  accounts  for  the  small  number  of  these 
broad  anomalies  which  were  available  for  interpretation  after 
numerical  filtering  was  applied. 

In  geophysics,  for  any  potential  field  anomaly 
measured  at  the  surface  there  are  basically  two  methods  by 
which  an  interpretation  for  the  possible  structure  causing 
the  anomaly  may  be  made.  The  word  possible  is  employed  here 
since  it  is  a  well-known  fact  that  potential  field  anomalies 
can  have  several  solutions,  all  equally  feasible. 

The  first  of  these  methods  is  that  which  employs 
downward  or  upward  continuation.  The  basic  concept  of  this 
method  is  that  one  calculates  a  set  of  coefficients  which 


are  then  convoluted  with  the  anomaly  data  such  that  the 
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anomaly  Is  continued  downward  to  the  depth  at  which  the  source 
is  located  or  an  anomaly  can  be  calculated  for  a  theoretical 
structure  at  a  certain  depth  and  then  continued  upward.  In 
the  former  method  the  calculated  anomaly  will  begin  to  oscil¬ 
late  once  the  approximate  location  of  the  true  depth  to  the 
source  is  passed.  In  the  latter  method  the  anomaly  continued 
upward  to  the  surface  must  match  the  observed  anomaly.,  then 
the  theoretical  structure  is  said  to  be  a  possible  source 
for  such  an  anomaly.  In  downward  continuation  one  can  project 
either  the  potential  or  the  field  downwards „  the  former 
provides  a  better  insight  as  to  the  shape  of  the  source. 

There  are  several  treatises  on  this  method  of  continuation ^ 
the  reader  is  referred  to  Peters  (19^9)  and  to  Dean  (1958) 
for  a  more  comprehensive  account  of  the  application  of  this 
method . 

The  second  method  is  that  of  geometric  construction 
of  a  cross  section  for  which  the  theoretical  anomaly  is  cal¬ 
culated.  The  shape  of  the  cross  section  as  well  as  its 
depth  are  then  altered  until  the  resulting  theoretical  anomaly 
approximates  as  accurately  as  possible  the  shape  of  the 
observed  anomaly. 

It  was  first  thought  that  the  interpretation  of 
the  seven  anomalies  found  by  numerical  filtering  might  best 
be  accomplished  by  calculating  the  potential  at  the  surface 
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using  the  filtered  vertical  component  of  the  field  and  down¬ 
wardly  continuing  this  potential  field  to  the  true  depth. 
However,  this  idea  was  abandoned  when  it  was  found,  as  pointed 
out  by  Dean  (1958),  that  the  data  errors  are  amplified  by  a 
factor  of  approximately  twice  the  difference  of  the  maximum 
adjacent  coefficients.  The  maximum  adjacent  coefficients  for 
the  depths  to  which  the  anomalies  were  to  be  projected  in 
this  analysis  differed  by  as  much  as  a  factor  of  lCr  in  some 
cases.  For  this  reason  the  anomalies  were  fitted  geometrically 
according  to  the  second  method  outlined  above. 

It  was  hypothesized  that  the  anomalies  were  due  to 
an  upwarping  or  downwarping  of  the  Curie  point  geotherm.  As 
shown  in  figure  13,  an  upwarping  of  the  geotherm  gives  rise 
to  a  negative  anomaly  in  the  earth’s  field  as  measured  at  the 
surface  and  a  downwarping  causes  a  positive  anomaly  to  appear 
at  the  surface  if  the  susceptibility  contrast  is  sufficient 
in  both  cases. 

The  seven  anomalies  were  then  fitted  geometrically 
with  theoretical  anomalies  using  the  above  facts  and  the 
well-known  theory1  for  the  calculation  of  the  vertical  mag¬ 
netic  attraction  of  a  distribution  of  magnetic  charge  induced 
by  the  earth’s  field  on  the  surface  of  a  two-dimensional 
structure.  The  following  assumptions  and  estimates  were  made 
in  employing  this  theory: 


This  theory,  as  outlined  by  Press  and  Ewing  (1952),  is  devel¬ 
oped  fully  in  appendix  3  for  the  two  types  of 
theoretical  structures  used  in  this  interpretation. 
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(i)  it  was  assumed  that  the  polarization  induced  by  the 
earth’s  field  on  the  body  was  in  the  vertical  direction  only 
since  the  majority  of  the  anomalies  were  relatively  symmetrical. 

(ii)  a  simplified  cross  section  symmetrical  about  its  central 
axis  was  assumed  for  each  structure,  again  because  of  the 
symmetry  of  the  observed  anomalies. 

(iii)  the  maximum  thickness  of  each  structure  was  assumed  to 
be  small  with  respect  to  its  width  (3  miles  as  compared  to 
100  miles  in  most  cases)  in  order  that  the  hypothesis  of  the 
source  being  a  warping  of  the  geotherm  be  plausible. 

(iv)  the  interpreter  must  estimate  the  lateral  extent  of 
each  anomaly  and  also  take  into  consideration  the  effects 
of  adjacent  anomalies  in  order  that  the  "end"  errors  be 
minimized . 

In  figures  14,  15,  16,  17,  18,  19,  and  20  are 
shown  the  final  two  dimensional  cross  sections  together 
with  the  corresponding  theoretical  and  observed  anomalies. 

These  interpretations  are  by  no  means  unique  solutions  for  the 
observed  anomalies;  but  they  are  thought  to  be  the  best  trial- 
and-error  solutions  under  the  assumptions  made  by  the  inter¬ 
preter  . 

In  order  to  investigate  the  possibility  that  the 
susceptibilities  of  continental  rocks  were,  in  general, 
greater  than  those  of  oceanic  rocks  as  suggested  by  the  inter¬ 
pretation  of  the  power  spectra,  the  susceptibility  contrasts 
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were  calculated  using  the  maximum  value  for  the  vertical  inten¬ 
sity  for  each  of  the  observed  and  theoretical  anomalies. 

The  values  for  the  vertical  component  of  the  earth’s  field 
were  taken  from  the  aeromagnetic  map  of  the  vertical  compon¬ 
ent  of  the  earth's  field  as  published  by  the  Hydrographic 
Office  of  the  United  States  Navy  (1955)*  Since  the  base  line, 
drawn  through  each  of  the  original  profiles  to  remove  the 
dipole  field,  was  not  absolutely  accurate  these  values  for 
susceptibility  contrast  can  only  be  taken  as  a  rough  estimate 
of  the  true  value  for  the  particular  structure  in  question. 

The  susceptibility  contrast  estimates  together  with  the  depth 
below  the  surface  of  the  earth  to  the  theoretical  structures 
are  summarized  in  table  3. 


Table  3 

Earth's  Anomaly  Features 


Region 

Flight  No. 

Vertical 
Component 
( oersteds ) 

No. 

Depthi 

(miles ) 

Susceptibility 
Contrast 
( c . g . s . )  x  108 

Oceanic 

10 

0.42 

A 

21.5 

5 

0.38 

B 

19.5 

10 

20 

0.48 

A 

18.5 

6 

23 

0.44 

A 

5.0 

15 

Continental  2 

O.58 

A 

22.5 

15 

0.58 

B 

25.5 

34 

7 

0.59 

A 

29.6 

15 

The  depths  for  continental  features  are  referred  to  the  land 

surface  while  those  for  oceanic  features  are  referred 
to  the  ocean  floor.  The  ocean  depths  are  only  approx¬ 
imate  and  were  obtained  from  the  "Advanced  Atlas  of 
Modern  Geography", (1962) . 
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The  calculations  of  susceptibility  contrast,  as 
summarised  in  table  3*  seem  to  support  the  idea  of  higher 
susceptibilities  for  continental  rocks  as  suggested  by  the 
interpretation  of  the  power  spectra  in  section  4.1.  In 
general,  the  depths  to  the  anomalous  structures,  which  have 
been  hypothesized  to  be  Curie  point  geotherm  structures  are 
slightly  larger  for  continental  features. 

The  anomaly  for  flight  23  does  not  agree  with  the 
pattern  set  forth  by  the  other  six  anomalies;  this  could  be 
due  to  the  fact  that  the  anomaly  location  is  quite  close  to 
the  terminal  point  of  the  flight,  Shemya  in  the  Aleutian 
Islands.  This  is  a  region  of  high  volcanic  activity  and 
the  temperature  gradient  here  may  be  much  higher  than  the 
mean  gradient  for  the  oceanic  crust,  thus  giving  rise  to  a 
shallow  Curie  point  geotherm.  The  anomaly  depth  is  also 
within  range  of  the  Mohorovicic  discontinuity  for  the 
oceanic  crust  (approximately  4  miles)  and  a  warping  of  the 
discontinuity  might  be  the  cause. 

Generally,  it  appears  that  the  depths  to  the  oceanic 
sources  are  well  below  the  average  depth  to  the  Mohorovicic 
discontinuity,  whereas,  those  of  the  continental  features 


68 


seem  to  be  in  the  discontinuity  depth  range  for  the  contin¬ 
ental  crust  (approximately  25  miles) .  Therefore,  it  is 
conceivable  that  the  Curie  point  geotherm  is  quite  closely 
associated  with  the  discontinuity  for  the  continental  crust; 
or,  perhaps  some  or  all  of  these  continental  anomalies  are 
due  to  anomalous  magnetic  features  of  the  discontinuity 
itself.  The  latter  interpretation  then  implies  that  the 
geotherm  depth  would  be  deeper  still  and  perhaps  the  suscepti¬ 
bility  contrast  there  is  insufficient  to  produce  an  anomaly 
on  the  earth’s  field  at  the  surface. 

The  accuracies  of  the  depths  and  susceptibility 
estimates  summarized  in  table  3  are  difficult  to  estimate. 
However,  it  is  safe  to  say  that  the  accuracy  of  the  depth 
estimates  is  better  than  that  of  the  susceptibility  contrasts. 
In  view  of  the  method  of  interpretation  and  the  method  of 
removing  the  long  period  dipole  field  trend  from  the  original 
data,  the  depth  estimates  are  probably  accurate  to  within 
±  3  miles  and  the  susceptibility  contrasts  to  within  ^  50  per 
cent  of  the  values  shown  in  table  3- 
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Locations  of  Anomalies  in  Figures  1,2,3 


69 . 


U 

>> 

0 

u 

> 

X? 

c. 

0 

> 

L. 

OJ 

£ 

V> 

0 

£_ 

03 

E 

0 

0 

o 

o 

o 

U) 

co 

n 

0 

c 

0 

„Q 

< 

JO 

< 

o 

I- 

o 


10 

cd 

it 


-o 


0 

> 

o 

_Q 

0 


CO 

Cu 

E 


-i-  > 


ANOMALY 'A' FOR  FLIGHT  2  FIGURE  14 
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ANOMALY'S'  FOR  FLIGHT  2  HGU 
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ANOMALY 'A'  FOR  FLIGHT  7  FIGURE  16 
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ANOMALY  A  FOR  FLIGHT  10  ~  FIGURE  17 
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ANOMALY 'B' FOR  FLIGHT  10  FIGURE  18 
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ANOMALY 'A'  FOR  FLIGHT  20  FIGURE  19 
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ANOMALY  X  FOR  FLIGHT  23  FIGURE  20 


V. 


CONCLUSION  AND  SUGGESTIONS  FOR  FURTHER  WORK 


In  view  of  the  fact  that  the  amount  of  data  avail¬ 
able  for  the  analyses  and  interpretations  involved  in  this 
thesis  was  small,  the  conclusions  made  here  must  be  termed 
suggestive  as  opposed  to  factual.  In  an  effort  to  analyze 
the  magnetic  character  of  the  earth’s  crust  as  a  whole,  one 
can  see  that  a  great  amount  of  data  would  be  required  in 
order  to  form  concrete  conclusions.  However,  the  results  of 
this  work  agree  favourably  with  those  of  the  majority  of 
research  in  this  field  and  one  can  say  that  the  interpreta¬ 
tions  set  forth  here  are  reasonably  accurate  and  entirely 
feasible  based  on  the  data  analyzed. 

As  pointed  out  in  section  4.1  there  is  good  agree¬ 
ment  between  the  results  of  the  power  spectral  analyses  and 
those  of  the  harmonic  analyses  of  Alldredge  and  Van  Voorhis 
(1961).  In  short,  both  indicate  a  non-magnetic  character 
for  the  greater  portion  of  the  mantle.  From  the  evidence  of 
the  weighted  spectra  for  continental  and  oceanic  regions  one 
might  say  that  spectral  analysis  might  prove  to  be  a  useful 
tool  in  differentiating  between  continental  and  oceanic  crusts 
in  areas  covered  by  ice,  for  example.  However,  the  evidence 
is  insufficient  to  indicate  a  definite  difference  in  the  amp¬ 
litudes  of  peaks  appearing  on  continental  and  oceanic  spectra 
further  work  along  these  lines  may  shed  some  light  on  this 
possibility . 
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In  general.,  the  Curie  point  depth  estimates  were 
deeper  for  the  continental  crust  than  those  found  by  Vacquier 
and  Affleck  (11  -  15  miles,  1941).  Moreover,  the  values 
found  in  the  work  of  this  thesis  were  in  the  range  of  20  -  25 
miles,  which  was  the  Curie  point  depth  range  predicted  for 
magnetite  by  Birch  (1955).  Birch’s  estimates  were  based  on 
calculations  using  temperature-depth  curves  for  various  models 
of  the  crust  (see  section  1.4).  In  order  to  check  the  hypothesis 
that  the  broad  anomalies  (periods  >  200  miles)  are  due  to 
anomalous  features  of  the  Curie  point  geotherm,  one  might 
investigate  areas  where  anomalous  heat  flows  have  been  recorded. 
In  this  way  one  could  determine  if  the  geotherm  depth  was 
greater  in  areas  of  lower  heat  flow  and  smaller  for  areas  of 
high  heat  flow. 

From  an  analysis  of  similar  magnetic  profiles  Serson 

_  O 

and  Hannaford  (1957)  predicted  susceptibilities  >  5  x  lO-3 

c.g.s.  for  the  oceanic  crustal  rocks  and  for  the  continental 
>  -2 

crustal  rocks,  /  5  x  10  c.g.s.  The  estimates  made  in  section 
4.2  agree  with  the  above  values  for  the  oceanic  crust  but  are 
smaller  by  a  factor  of  2  or  3  for  the  continental  crust.  The 
only  depth,  which  might  be  interpreted  as  a  Curie  point  geo¬ 
therm  depth,  given  by  Serson  and  Hannaford, was  that  for  the 
lower  boundary  of  the  layer  of  dipoles  representing  the 
oceanic  magnetic  crust.  The  depth  given  was  ^  30  km  (approx¬ 
imately  20  miles),  which  is  in  general  agreement  with  the 
oceanic  values  of  section  4.2.  Along  these  channels,  it 
might  prove  fruitful  to  carry  out  more  accurate  surveys 


, 
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at  altitudes  cf  8000  -  10,000  feet  in  order  to  provide  more 
reliable  estimates  of  the  susceptibility  contrasts  at 
greater  depths  in  the  crust.  In  this  way,  one  might  obtain 
a  better  knowledge  of  the  magnetic  character  of  the  rocks 
at  these  depths  and  hence  a  better  understanding  of  the 
general  nature  of  the  rocks  themselves. 

In  future,  the  following  suggestions  may  prove 
helpful  in  planning  surveys  for  long  aeromagnetic  profiles 
from  which  broad  anomalies  are  to  be  extracted  and  interpreted: 

(i)  maintain  as  high  an  accuracy  as  possible  with  respect 

to  aircraft  speed  and  navigation  in  order  that  anomalies 
can  be  accurately  checked  with  surface  and  subsurface 
geology, 

(ii)  survey  those  areas  in  which  a  good  control  for  the 
dipole  field  trend  can  be  obtained,  if  possible,  to 
insure  maximum  accuracy  for  anomaly  amplitudes, 

(iii)  survey  some  of  the  chosen  areas  more  than  once  to  see 

if  magnetic  fluctuations  affect  the  records  drastically, 

(iv)  fly  global  flights1  (great  circles),  if  possible, 
avoiding  sharp  turns,  in  order  to  investigate  the 
magnetic  properties  of  the  core-mantle  boundary. 


1 


Such  flights  are  being  made  at  present  by  project  Magnet  of 
the  U.  S.  Hydrographic  Office. 
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APPENDIX  I 


The  Spectral  Window 


In  order  to  recognize  the  need  for  the  application 
of  a  spectral  window  in  order  to  obtain  a  power  spectrum 
for  a  set  of  discrete  data,  we  shall  first  consider  the  con¬ 
tinuous  series  rather  than  the  discrete  series. 

Given  a  continuous  record  of  finite  length,  it  is 
clear  that  we  cannot  estimate  the  autocovariance  function 
C(A)  for  arbitrarily  long  lags.  Indeed,  no  estimate  can  be 
made  for  lags  longer  than  the  record.  Thus,  in  place  of 


C(A) 


lim 

X-*t» 


1  X/2 

X  -X/2 


H(x) 


H(x+A)  dx 


for  all  values  of  A,  we  will  have 

i  ( Xn  -  I  A  I  )  /2 

cnn ( A)  =  — - -  f  H(x-A/2) -H(x+A/2)dx 

VIM  -(xn-|A|)/2 

=  coc/~^ 

only  for  |  A  |  ^  X  <  X  ,  where  Xn  is  the  length  of  the 

record,  and  Xm  is  the  maximum  lag  to  be  used.  CQq(A)  will 
be  called  the  "apparent  autocovariance  function", 
average  value  is  C(A)  for  |  A  |  ^  X  . 


since  its 
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For  the  continuous  series,  the  power  spectral 
estimates  are  usually  derived  from  a  "modified  apparent 
autocovariance  function"  by  Fourier  transformation. 

While  the  modified  apparent  autocovariance  func¬ 
tions,  which  are  obtained  by  multiplying  the  apparent  auto¬ 
covariance  function  by  suitable  even  functions  of  A,  are 
often  far  from  being  respectable  estimates  of  the  true 
autocovariance  function,  their  transforms  are  very  respect¬ 
able  estimates  of  "smoothed"  values  of  the  true  spectral 
density. 

Let  D^(A)  be  a  prescribed  even  function  of  A, 
subject  to  the  restriction  D^(0)  =  1,  and  Ib(A)  =  0  for 
|  A  |  >  X  ,  (where  i  =  0,  1,  2,  3*  depending  upon  the  shape 

of  D^(A)  for  |  A  [  <  Xm) ,  and  let  the  corresponding  modified 

apparent  autocovariance  function  be  defined  by 

C.(A)  =  D.(A)  *  C00(A) 

D.(A)  may  be  regarded  as  a  window  of  variable  transmission 
which  modifies  the  values  of  Cqq(A)  differently  for  differ¬ 
ent  lags.  It  is  therefore  natural  to  call  D^(A)  a  "lag 
window . " 

It  is  clear  that  Ch(A)  =  0  for  |  A  |  >  Xm 

although  C00( A)  was  not  defined  there.  Because  Ch(A)  is 
defined  for  all  values  of  A,  it  has  a  perfectly  definite 
Fourier  transform  P^(f),  which  should  satisfy  the  relation. 


i 
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P±(f)  -  Q±(f)  *  P00(f) 
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where  Q^(f)  is  the  Fourier  transform  of  D^(A),  the  asterisk 
indicates  convolution,  and  P  ^(f)  is  the  Fourier  transform 
of  Cqq(^) *  However,  P  ^(f)  is  not  determinate  because  CQq(A) 
is  not  specified  for  |A  |  >  Xm  (and  its  definition  cannot  be 
directly  extended  beyond  |  A  |  =  3^).  Nevertheless,  since 

ave  [C±( A) ]  =  D.(A)  •  C( A) 

where  C(A)  is  the  true  autocovariance  function,  it  follows 
that 

ave  [P±(f ) ]  =  Q.(f)  *  P(f) 

where  P(f)  is  the  true  power  spectrum,  that  is,  the  Fourier 
transform  of  C(A).  The  average  may  be  thought  of  as  along 
space.  This  corresponds  to  replacing  H(x)  by  H(x-A) ,  thus 
changing  the  stretch  of  H(x)  which  is  observed,  and  then 
averaging  over  A.  The  corresponding  relation, 

00 

ave  [P.(f1)]  =  f"  Q±(f1-f)  •  P(f)  df 

exhibits  the  average  value  of  P^(f^)  as  a  smoothing  (average- 
over-frequency)  of  the  true  power  spectrum  density  P(f) 
over  frequencies  "near"  f^,  with  weights  proportional  to 
Q^(f-^-f).  Then  P^(f  )  may  be  said  to  be  the  collected 
impression  of  the  true  power  spectrum  P(f)  obtained  through 


A5 


a  window  of  variable  transmission  Q^(f^-f).  It  is  therefore 
natural  to  call  Q1(f)  the  "spectral  window"  corresponding  to 
the  "lag  window"  D^(A). 

In  the  case  of  discrete  series  we  again  run  into 
the  problem  of  being  unable  to  apply  arbitrarily  long  lags. 

As  a  consequence  we  have  to  deal  with  finite  Fourier  series 
transformations  rather  than  with  infinite  Fourier  integral 
transformations,  however  the  correspondence  between  multi¬ 
plication  and  convolution  persists. 

The  application  of  a  spectral  window  in  the  discrete 
series  case  is  slightly  different.  In  the  continuous  case 
the  window  is  introduced  as  a  lag  window  before  cosine  trans¬ 
formation.,  in  order  to  avoid  transformation  after  convolution. 
In  the  discrete  case,,  since  convolution  is  not  difficulty 
convolution  after  transformation  is  convenient,  and  the 
window  is  shaped  after  cosine  transformation. 

The  reason  for  this  change  of  procedure  is  that 
the  D^(A),  for  i  >  1,  are  finite  sums  of  cosines,  so  that 
their  transforms  are  simply  sums  of  spikes  (Dirac  delta- 
functions)  at  the  appropriate  spacing.  Hence  convolution 
becomes  only  smoothing  with  a  set  of  weights  whose  values 
depend  on  the  shape  of  the  spectral  window  desired. 

In  choosing  the  proper  lag  window,  one  would  like 
to  concentrate  the  main  lobe  of  Qj_(f)  near  f  =  0,  keeping 
the  side  lobes  as  low  as  feasible.  In  order  to  concentrate 
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the  main  lobe  one  must  make  D.-l(a),  flat  and  rather  blocky. 

In  order  to  reduce  the  side  lobes,  however,  one  has  to  make 
D^(A)  smooth  and  gently  changing.  Since  D^(A)  must  vanish 
for  |  A  |  >  Xm  a  compromise  must  be  made.  Such  a  compromise 

is  described  below. 


f c|  u  r  e  Z 

The  lag  window  shown  in  figure  1  has  the  form 
D  (X)  =  1/2(1  +  cos  |X|  <  Xm 

=  0  I x  I  >  xm 

and  the  corresponding  spectral  window,  shown  in  figure  2, 
is  written 

^(f)  =  1/2  Qoff)  +  l/4[0o(f+  ^-)  +  Q0(f-  2^)] 

Hence  the  convolution  with  this  particular  spectral  window 
means  smoothing  with  weights  0.25,  0.5*  0.25.  This  form  of 
smoothing  is  known  as  "hanning",  after  the  Austrian  meteorologist 


Julius  von  Hann. 
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APPENDIX  II 

Numerical  Calculation  of  the  Convolution  Integral 

The  convolution  of  two  functions  Is  accomplished 
by  means  of  the  convolution  Integral  given  as, 

Y(6)  =  f°°  W(x)  X ( 6-x )  dx  .  (1) 

CO 

For  discrete  data  the  above  integral  degenerates 

to  a  sum, 

N 

Y(6)  =  2  W(a)  X(  6-a)  Aa  .  (2) 

a=l 

where  Y  -  output 

X  -  input 

W  -  weighting  function 
6  -  lag 

Ac  -  sampling  interval 

The  equation  (2)  is  easily  programmed  for  the  digital 
computer  since  it  may  be  written  as, 

=  W-^X-^Aa 

Y2  =  (WXX2  +  W2X]_)Aa 

YN  =  (W1XN  +  W2XN-l+,,,+  WN-1X2  +  WNXl^Aa 

that  is,  the  weighting  moves  from  left  to  right  on  the  profile, 
moving  one  sampling  interval  for  each  computation  of  one  out¬ 
put  point. 


APPENDIX  III 


Calculation  of  Theoretical  Anomalies  Due  to  Warping  of  the 

Curie  Point  Oeotherm 

These  calculations  were  made  assuming 

(i)  a  distribution  of  magnetic  charge  induced  by  the  earth 3 s 
field  on  the  surface  of  a  two-dimensional  structure,, 

(ii)  polarization  of  the  structure  in  the  vertical  direction 
only. 


Figure  1 

Prom  figure  1,  the  potential  at  the  point  P ,  due 
to  the  line  element  of  magnetic  charge  ds,  is  given  by 

2  k  Zq  cos  0  loger  ds 

where  Zq  is  the  vertical  component  of  the  earth’s  field,  $  is 
the  angle  of  inclination  of  the  line  element  ds,  and  k  is  the 
magnetic  susceptibility  contrast  for  the  structure. 
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The  field  at  P  due  to  the  line  element  is 
(2  k  ZQ  cos  0  ds)/r 

and  the  vertical  component  of  this  field  is 

2  k  ZQ  cos  0  ds 

where  r  and  9  are  the  magnitude  and  direction  of  the  radius 
vector  from  the  point  of  observation,  P,  to  the  line  element  ds . 

To  obtain  the  expression  for  the  field  due  to  the 
entire  line  distribution  with  inclination  0  we  must  integrate 
over  the  contour  C,  giving, 

AZ  =  2  k  Jc  (ZQ  cos  0  -5iS-2)  ds 

From  figure  1, 

ds  sin  9  =  rd@  cos  0  +  dr  sin  0 

.  .  AZ  =  2  k  fg  ZQ  cos  0  [sin  0(&r/v)  -  cos  0  d0]  . ..(l) 

To  obtain  the  field  component  for  the  entire 
structure  the  contour  integral  must  be  carried  out  over  the 
entire  cross  section  of  the  structure.  If  we  represent  the 
cross  section  by  a  sequence  of  line  segments,  one  of  which 
is  shown  in  figure  1,  (l)  can  be  integrated  over  each  seg¬ 
ment,  and  the  results  can  be  summed  over  all  segments  to  give 

=  2  k  2  Z0  cos  pfJcoB  pri(01  -  01+1)  +  sin  loge£|±i)  ]  .  .  .  (2) 


AZ 


.  .  u 

:'<0  i: -J J..  ■  ■  9  ;ii.  i-y  ,_  '  HO 


A10 


4-  U 

where  is  the  inclination  of  the  i  segment,  and  r^,  0^. 

r . , , ,  0,  ,  _  define  the  radius  vectors  from  the  point  of 
l+l  l+l 

4“  V~\ 

observation  to  the  end  points  of  the  i  segment. 

A  downwarping  of  the  Curie  point  geotherm  may  be 
represented  by  the  following  two-dimensional  structure: 


e-9~a 

r  3  4 


e*©-©. 

8  2.  5 


From  equation  2  we  have  that 


AZ  =  AZ23  +  +  AZ^  +  AZ^2 


and  AZ 


?  ro 

23  =  2  k  ZQ[  cos  $(e2-Q^)  +  cos  0  sin  pf  loge 


^^34  -  ^  k  Z0( 63-0^) 

AZ45  =  2  k  ZQ[cos2  ^(0^-0^)  -  cos  0  sin  0  loge  ^  ] 

AZ52  =  2  k  Zq(05-02) 

Collection  of  like  terms  gives 

AZ  =  2  k  ZQ[sin2  ^(0T-0B)  +  cos  $  sin  $  loge  ] 


■ 
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where  0T  =  angle  suspended  by  the  top  horizontal  line  segment 
=  angle  suspended  by  the  bottom  horizontal  line 
segment . 

Similarly,  one  may  represent  an  upwarplng  of  the 
geotherm  by  the  following  two-dimensional  structure: 


Upon  applying  the  mathematical  treatment  used  In 
the  previous  case,  one  obtains  the  following  expression 
for  the  vertical  component  of  the  field  due  to  the  theoretical 
structure: 


=  2k  ZQ[sin2  ^(0T-@B)  +  sin  $  cos  0  log 
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